##Install Packages if Needed
if (!require("ggplot2")) install.packages("ggplot2")
if (!require("cowplot")) install.packages("cowplot")
if (!require("lubridate")) install.packages("lubridate")
Loading required package: lubridate
Attaching package: ‘lubridate’
The following object is masked from ‘package:cowplot’:
stamp
The following objects are masked from ‘package:base’:
date, intersect, setdiff, union
if (!require("corrplot")) install.packages("corrplot")
if (!require("lme4")) install.packages("lme4")
if (!require("vegan")) install.packages("vegan")
if (!require("DHARMa")) install.packages("DHARMa")
if (!require("effectsize")) install.packages("effectsize")
if (!require("emmeans")) install.packages("emmeans")
if (!require("lmerTest")) install.packages("lmerTest")
if (!require("tidyr")) install.packages("tidyr")
if (!require("dplyr")) install.packages("dplyr")
if (!require("plotrix")) install.packages("plotrix")
Loading required package: plotrix
Warning: package ‘plotrix’ was built under R version 4.3.2
if (!require("Rmisc")) install.packages("Rmisc")
if (!require("ggpubr")) install.packages("ggpubr")
Loading required package: ggpubr
Attaching package: ‘ggpubr’
The following object is masked from ‘package:plyr’:
mutate
The following object is masked from ‘package:cowplot’:
get_legend
##Load Packages
library(ggplot2) #Required for ggplots
library(cowplot) #Required for plotting panel figures
library(lubridate) #Required for date and time formatting
library(corrplot) #Required for correlation plot
library(lme4) #Required for mixed effects modeling
library(vegan) #Required for multivariate analysis PERM
library(DHARMa) #Required to check residuals of mixed effects models
library(effectsize) #Required for eta_squared effect sizes
library(emmeans) #Required for pairwise comparisons
library(lmerTest) #Required for p values with lme4 model summaries
library(tidyr) #Required for pivot_wider function
library(dplyr) #Required for group_by function
library(plotrix) #Required for Standard error function
library(Rmisc) #Required for SummarySE function
library(ggpubr) #Required for stat brackets
#Note: Run "Graphing Parameters" section from 01_ExperimentalSetup.Rmd file
##Growth
Length<-read.csv("Data/LengthFull.csv", header=TRUE)
Dates<-read.csv("Data/BonaireDates.csv", header=TRUE)
##Set date variables
Dates$Date<-as.Date(Dates$Date, format='%m/%d/%Y')
##Set factor variables
Length$Site<-factor(Length$Site, levels=c("KL", "SS"))
Length$Genotype<-factor(Length$Genotype, levels=c("AC8", "AC10", "AC12"))
Length$Orig<-factor(Length$Orig, levels=c("N", "T"))
Length$Origin<-factor(Length$Origin, levels=c("Native", "Transplant"))
##Add a Sample Set Variable
Length$Set<-paste(Length$Site, Length$Genotype, Length$Orig, sep=".")
Length$Set<-factor(Length$Set)
##Add Site.Orig variable
Length$Site.Orig<-paste(Length$Site, Length$Orig, sep=".")
Length$Site.Orig<-factor(Length$Site.Orig, levels=c("KL.N", "KL.T", "SS.N", "SS.T"))
##Bleaching Metrics
#Note: Physiological metrics calculated in 02_PhysiologyMetrics.R file
Thermal<-read.csv("Outputs/ThermData.csv", header=TRUE)
##Set factor variables
Thermal$TimeP<-factor(Thermal$TimeP, levels=c("IN", "TP1", "TP2", "TP3", "TP4"))
Thermal$Site<-factor(Thermal$Site, levels=c("KL", "SS"))
Thermal$Genotype<-factor(Thermal$Genotype, levels=c("AC8", "AC10", "AC12"))
Thermal$Orig<-factor(Thermal$Orig, levels=c("N", "T"))
Thermal$Origin<-factor(Thermal$Origin, levels=c("Native", "Transplant"))
##Subset Dates for TP 1 and TP 2 Growth
Growth.Dates_T1.2<-subset(Dates, Task=="Growth" & TimeP=="TP1" |
Task=="Growth" & TimeP=="TP2" )
##Convert to Wide format
Growth.Dates_T1.2<-Growth.Dates_T1.2 %>% pivot_wider(names_from = TimeP, values_from = Date)
##Calculate Number of Days
#TP 1 to TP2
Growth.Dates_T1.2$TP1.2_days<-time_length(Growth.Dates_T1.2$TP2-Growth.Dates_T1.2$TP1, unit="days")
##Merge Length Data with Growth Dates
#Merges by Site column
#Adds TP1.2_days column
Growth_T1.2<-merge(Length, Growth.Dates_T1.2[,c(-2)], all.x=TRUE)
##Calculate Linear Extension (cm)
Growth_T1.2$Ext_cm<-Growth_T1.2$TL.TP2_cm-Growth_T1.2$TL.TP1_cm
##Calculate Growth_T1.2 Rate (cm/day)
Growth_T1.2$TLE_cm.day<-Growth_T1.2$Ext_cm/Growth_T1.2$TP1.2_days
##Remove Corals not Measured in T1 to T2
Growth_T1.2<-Growth_T1.2 %>% drop_na(TLE_cm.day)
##Sample Size per Set
Growth_T1.2 %>%
dplyr::group_by(Set) %>%
dplyr::summarize(SampleSize = length(Set))
##Check for Outliers
ggplot(Growth_T1.2, aes(x=Set, y=TLE_cm.day)) +
geom_boxplot(alpha=0.5, shape=2, outlier.shape = NA)+
geom_jitter(shape=16, position=position_jitter(0.1))+
theme(axis.text.x = element_text(angle = 90))
Potential outliers have been re-measured and therefore will be retained. n=10 per Site, Genotype, Origin group
#Check normality
hist(Growth_T1.2$TLE_cm.day)
shapiro.test(Growth_T1.2$TLE_cm.day)
Shapiro-Wilk normality test
data: Growth_T1.2$TLE_cm.day
W = 0.92803, p-value = 7.262e-06
#Not normal
hist(log(Growth_T1.2$TLE_cm.day+1))
shapiro.test(log(Growth_T1.2$TLE_cm.day+1))
Shapiro-Wilk normality test
data: log(Growth_T1.2$TLE_cm.day + 1)
W = 0.94335, p-value = 7.248e-05
##Still not normal
##Compare generalized linear mixed effects models
TLE_T1.2.glmr.gaus<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T1.2, family=gaussian(link="identity"))
Warning: calling glmer() with family=gaussian (identity link) as a shortcut to lmer() is deprecated; please call lmer() directly
TLE_T1.2.glmr.gaus.l<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T1.2, family=gaussian(link="log"))
TLE_T1.2.glmr.gaus.i<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T1.2, family=gaussian(link="inverse"))
TLE_T1.2.glmr.gam.l<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T1.2, family=Gamma(link="log"))
TLE_T1.2.glmr.gam.i<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T1.2, family=Gamma(link="inverse"))
TLE_T1.2.glmr.inv<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T1.2, family=inverse.gaussian(link="1/mu^2"))
boundary (singular) fit: see help('isSingular')
TLE_T1.2.glmr.inv.l<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T1.2, family=inverse.gaussian(link="log"))
TLE_T1.2.glmr.inv.i<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T1.2, family=inverse.gaussian(link="inverse"))
AIC(TLE_T1.2.glmr.gaus, TLE_T1.2.glmr.gaus.l, TLE_T1.2.glmr.gaus.i, TLE_T1.2.glmr.gam.l, TLE_T1.2.glmr.gam.i, TLE_T1.2.glmr.inv, TLE_T1.2.glmr.inv.l, TLE_T1.2.glmr.inv.i)
Gamma distribution with log-link has the lowest AIC.
##Check residuals
TLE_T1.2.glmr.gam.l_res <- simulateResiduals(fittedModel = TLE_T1.2.glmr.gam.l, plot = F)
plot(TLE_T1.2.glmr.gam.l_res)
Compare residuals across models
##Check residuals
TLE_T1.2.glmr.gaus.l_res <- simulateResiduals(fittedModel = TLE_T1.2.glmr.gaus.l, plot = F)
plot(TLE_T1.2.glmr.gaus.l_res)
Compare with log+1 transformed lmer model
##Model with log +1 transformation
#Function of Site and Origin, with Genotype as a Random effect
#Interaction between Site and Origin
TLE_T1.2.lme<-lmer(log(TLE_cm.day+1)~Origin*Site+(1|Genotype), data=Growth_T1.2)
##Check residuals
TLE_T1.2.lme_res <- simulateResiduals(fittedModel = TLE_T1.2.lme, plot = F)
plot(TLE_T1.2.lme_res)
Better residuals with the log+1 transformed lmer model.
##Model results
summary(TLE_T1.2.lme)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(TLE_cm.day + 1) ~ Origin * Site + (1 | Genotype)
Data: Growth_T1.2
REML criterion at convergence: -233.2
Scaled residuals:
Min 1Q Median 3Q Max
-2.0026 -0.5991 -0.1118 0.5959 2.4282
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.002965 0.05445
Residual 0.006632 0.08144
Number of obs: 120, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.125179 0.034778 2.684719 3.599 0.04400 *
OriginTransplant -0.007538 0.021027 114.000000 -0.359 0.72063
SiteSS 0.062096 0.021027 114.000000 2.953 0.00382 **
OriginTransplant:SiteSS -0.035167 0.029737 114.000000 -1.183 0.23943
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) OrgnTr SiteSS
OrgnTrnspln -0.302
SiteSS -0.302 0.500
OrgnTrn:SSS 0.214 -0.707 -0.707
eta_squared(TLE_T1.2.lme)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-------------------------------------------
Origin | 0.02 | [0.00, 1.00]
Site | 0.07 | [0.01, 1.00]
Origin:Site | 0.01 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
TLE_T1.2.res<-data.frame(summary(TLE_T1.2.lme)$coefficients[-1,])
TLE_T1.2.res$Predictor<-c("Origin", "Site", "Origin x Site")
TLE_T1.2.res$EtaSq<-c(eta_squared(TLE_T1.2.lme)$Eta2)
TLE_T1.2.res$Response<-rep("Growth", nrow(TLE_T1.2.res))
TLE_T1.2.res$Timepoint<-rep("TP1_2", nrow(TLE_T1.2.res))
Effect size of Origin for each Site
##KL
TLE_T1.2.lme_KL<-lmer(log(TLE_cm.day+1)~Origin+(1|Genotype), data=Growth_T1.2[which(Growth_T1.2$Site=="KL"),])
summary(TLE_T1.2.lme_KL)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(TLE_cm.day + 1) ~ Origin + (1 | Genotype)
Data: Growth_T1.2[which(Growth_T1.2$Site == "KL"), ]
REML criterion at convergence: -115.1
Scaled residuals:
Min 1Q Median 3Q Max
-1.7586 -0.6776 -0.1453 0.5574 2.2922
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.001699 0.04122
Residual 0.006724 0.08200
Number of obs: 60, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.125179 0.028116 2.712715 4.452 0.026 *
OriginTransplant -0.007538 0.021173 56.000006 -0.356 0.723
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.377
eta_squared(TLE_T1.2.lme_KL)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 2.26e-03 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##SS
TLE_T1.2.lme_SS<-lmer(log(TLE_cm.day+1)~Origin+(1|Genotype), data=Growth_T1.2[which(Growth_T1.2$Site=="SS"),])
summary(TLE_T1.2.lme_SS)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(TLE_cm.day + 1) ~ Origin + (1 | Genotype)
Data: Growth_T1.2[which(Growth_T1.2$Site == "SS"), ]
REML criterion at convergence: -114.6
Scaled residuals:
Min 1Q Median 3Q Max
-2.05792 -0.57732 -0.09493 0.72022 2.38351
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.004144 0.06437
Residual 0.006600 0.08124
Number of obs: 60, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.18727 0.04002 2.30548 4.680 0.0323 *
OriginTransplant -0.04271 0.02098 56.00000 -2.036 0.0465 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.262
eta_squared(TLE_T1.2.lme_SS)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.07 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save results
TLE_T1.2.site.res<-data.frame(rbind(summary(TLE_T1.2.lme_KL)$coefficients[-1,],
summary(TLE_T1.2.lme_SS)$coefficients[-1,]))
TLE_T1.2.site.res$Predictor<-c("KL Origin", "SS Origin")
TLE_T1.2.site.res$EtaSq<-c(eta_squared(TLE_T1.2.lme_KL)$Eta2, eta_squared(TLE_T1.2.lme_SS)$Eta2)
TLE_T1.2.site.res$Response<-rep("Growth", nrow(TLE_T1.2.site.res))
TLE_T1.2.site.res$Timepoint<-rep("TP1_2", nrow(TLE_T1.2.site.res))
##Combine results
TLE_T1.2.res<-rbind(TLE_T1.2.res, TLE_T1.2.site.res)
##Summary statistics by Site and Origin
TP1.2_Growth.sum<-summarySE(Growth_T1.2, measurevar="TLE_cm.day", groupvars=c("Site", "Origin", "Site.Orig"), na.rm=TRUE)
##Plot Average Growth Rate across Treatments
TP1.2_Growth.plot<-ggplot(TP1.2_Growth.sum, aes(x=Site, y=TLE_cm.day, colour=Site.Orig)) +
scale_colour_manual(values=Orig.colors.o)+
geom_errorbar(aes(ymin=TLE_cm.day-se, ymax=TLE_cm.day+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Origin", y=expression(paste('Growth Rate (cm day'^-1*")")), colour="Origin")+
ylim(0, 0.25)+
annotate("text", x=2, y=0.24, label="*", size=sig.sz, fontface="bold")+
geom_bracket(xmin=1, xmax=2, y.position=0.01, label.size=levels.sz, label="**", inherit.aes = FALSE); TP1.2_Growth.plot
##Summary statistics by Genotype and Origin
TP1.2_Growth.sum.SS<-summarySE(subset(Growth_T1.2, Site=="SS"), measurevar="TLE_cm.day", groupvars=c("Genotype", "Origin"), na.rm=TRUE)
TP1.2_Growth.sum.SS
##Calculate Percent Difference between Native and Transplant by Genotype
TP1.2_Growth.SS<-TP1.2_Growth.sum.SS[which(TP1.2_Growth.sum.SS$Origin=="Native"), c(1,4)]
names(TP1.2_Growth.SS)[2]<-paste(names(TP1.2_Growth.SS)[2], "N", sep="_")
TP1.2_Growth.SS$TLE_cm.day_T<-TP1.2_Growth.sum.SS$TLE_cm.day[which(TP1.2_Growth.sum.SS$Origin=="Transplant")]
TP1.2_Growth.SS$Perc.Inc<-((TP1.2_Growth.SS$TLE_cm.day_N-TP1.2_Growth.SS$TLE_cm.day_T)/TP1.2_Growth.SS$TLE_cm.day_T)*100
TP1.2_Growth.SS
##Mean and Standard Error
mean(TP1.2_Growth.SS$Perc.Inc)
[1] 47.19841
std.error(TP1.2_Growth.SS$Perc.Inc)
[1] 25.89777
##Subset Dates for TP 3 and TP 4 Growth
Growth.Dates_T3.4<-subset(Dates, Task=="Growth" & TimeP=="TP3" |
Task=="Growth" & TimeP=="TP4" )
##Convert to Wide format
Growth.Dates_T3.4<-Growth.Dates_T3.4 %>% pivot_wider(names_from = TimeP, values_from = Date)
##Calculate Number of Days
#TP 3 to TP4
Growth.Dates_T3.4$TP3.4_days<-time_length(Growth.Dates_T3.4$TP4-Growth.Dates_T3.4$TP3, unit="days")
##Merge Length Data with Growth Dates
#Merges by Site column
#Adds TP3.4_days column
Growth_T3.4<-merge(Length, Growth.Dates_T3.4[,c(-2)], all.x=TRUE)
##Calculate Linear Extension (cm)
Growth_T3.4$Ext_cm<-Growth_T3.4$TL.TP4_cm-Growth_T3.4$TL.TP3_cm
##Calculate Growth_T3.4 Rate (cm/day)
Growth_T3.4$TLE_cm.day<-Growth_T3.4$Ext_cm/Growth_T3.4$TP3.4_days
##Remove Corals not Measured in T1 to T2
Growth_T3.4<-Growth_T3.4 %>% drop_na(TLE_cm.day)
##Sample Size per Set
Growth_T3.4 %>%
dplyr::group_by(Set) %>%
dplyr::summarize(SampleSize = length(Set))
##Check for Outliers
ggplot(Growth_T3.4, aes(x=Set, y=TLE_cm.day)) +
geom_boxplot(alpha=0.5, shape=2, outlier.shape = NA)+
geom_jitter(shape=16, position=position_jitter(0.1))+
theme(axis.text.x = element_text(angle = 90))
#Check normality
hist(Growth_T3.4$TLE_cm.day)
shapiro.test(Growth_T3.4$TLE_cm.day)
Shapiro-Wilk normality test
data: Growth_T3.4$TLE_cm.day
W = 0.92767, p-value = 6.904e-06
#Not normal
hist(log(Growth_T3.4$TLE_cm.day+1))
shapiro.test(log(Growth_T3.4$TLE_cm.day+1))
Shapiro-Wilk normality test
data: log(Growth_T3.4$TLE_cm.day + 1)
W = 0.94719, p-value = 0.0001351
##Still not normal
##Compare generalized linear mixed effects models
TLE_T3.4.glmr.gaus<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T3.4, family=gaussian(link="identity"))
Warning: calling glmer() with family=gaussian (identity link) as a shortcut to lmer() is deprecated; please call lmer() directly
TLE_T3.4.glmr.gaus.l<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T3.4, family=gaussian(link="log"))
TLE_T3.4.glmr.gaus.i<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T3.4, family=gaussian(link="inverse"))
TLE_T3.4.glmr.gam.l<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T3.4, family=Gamma(link="log"))
TLE_T3.4.glmr.gam.i<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T3.4, family=Gamma(link="inverse"))
TLE_T3.4.glmr.inv<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T3.4, family=inverse.gaussian(link="1/mu^2"))
TLE_T3.4.glmr.inv.l<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T3.4, family=inverse.gaussian(link="log"))
TLE_T3.4.glmr.inv.i<-glmer(TLE_cm.day~Origin*Site+(1|Genotype), data=Growth_T3.4, family=inverse.gaussian(link="inverse"))
AIC(TLE_T3.4.glmr.gaus, TLE_T3.4.glmr.gaus.l, TLE_T3.4.glmr.gaus.i, TLE_T3.4.glmr.gam.l, TLE_T3.4.glmr.gam.i, TLE_T3.4.glmr.inv, TLE_T3.4.glmr.inv.l, TLE_T3.4.glmr.inv.i)
Gamma distribution with log-link has the lowest AIC.
##Check residuals
TLE_T3.4.glmr.gam.l_res <- simulateResiduals(fittedModel = TLE_T3.4.glmr.gam.l, plot = F)
plot(TLE_T3.4.glmr.gam.l_res)
Compare residuals across models
##Check residuals
TLE_T3.4.glmr.gaus.l_res <- simulateResiduals(fittedModel = TLE_T3.4.glmr.gaus.l, plot = F)
plot(TLE_T3.4.glmr.gaus.l_res)
Compare with log+1 transformed lmer model
##Model with log +1 transformation
#Function of Site and Origin, with Genotype as a Random effect
#Interaction between Site and Origin
TLE_T3.4.lme<-lmer(log(TLE_cm.day+1)~Origin*Site+(1|Genotype), data=Growth_T3.4)
##Check residuals
TLE_T3.4.lme_res <- simulateResiduals(fittedModel = TLE_T3.4.lme, plot = F)
plot(TLE_T3.4.lme_res)
Better residuals with the log+1 transformed lmer model.
##Model results
summary(TLE_T3.4.lme)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(TLE_cm.day + 1) ~ Origin * Site + (1 | Genotype)
Data: Growth_T3.4
REML criterion at convergence: -341.9
Scaled residuals:
Min 1Q Median 3Q Max
-2.27432 -0.60263 -0.01933 0.56722 2.45603
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.002008 0.04481
Residual 0.002574 0.05074
Number of obs: 120, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.16891 0.02748 2.38972 6.147 0.016262 *
OriginTransplant -0.02845 0.01310 114.00000 -2.172 0.031936 *
SiteSS -0.05868 0.01310 114.00000 -4.479 1.79e-05 ***
OriginTransplant:SiteSS 0.06284 0.01853 114.00000 3.392 0.000954 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) OrgnTr SiteSS
OrgnTrnspln -0.238
SiteSS -0.238 0.500
OrgnTrn:SSS 0.169 -0.707 -0.707
eta_squared(TLE_T3.4.lme)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-------------------------------------------
Origin | 9.01e-04 | [0.00, 1.00]
Site | 0.07 | [0.01, 1.00]
Origin:Site | 0.09 | [0.02, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
TLE_T3.4.res<-data.frame(summary(TLE_T3.4.lme)$coefficients[-1,])
TLE_T3.4.res$Predictor<-c("Origin", "Site", "Origin x Site")
TLE_T3.4.res$EtaSq<-c(eta_squared(TLE_T3.4.lme)$Eta2)
TLE_T3.4.res$Response<-rep("Growth", nrow(TLE_T3.4.res))
TLE_T3.4.res$Timepoint<-rep("TP3_4", nrow(TLE_T3.4.res))
Effect size of Origin for each Site
##KL
TLE_T3.4.lme_KL<-lmer(log(TLE_cm.day+1)~Origin+(1|Genotype), data=Growth_T3.4[which(Growth_T3.4$Site=="KL"),])
summary(TLE_T3.4.lme_KL)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(TLE_cm.day + 1) ~ Origin + (1 | Genotype)
Data: Growth_T3.4[which(Growth_T3.4$Site == "KL"), ]
REML criterion at convergence: -155.5
Scaled residuals:
Min 1Q Median 3Q Max
-2.13557 -0.60999 -0.09472 0.50939 2.13461
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.002615 0.05113
Residual 0.003233 0.05686
Number of obs: 60, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.16891 0.03129 2.23939 5.398 0.0254 *
OriginTransplant -0.02845 0.01468 56.00000 -1.938 0.0577 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.235
eta_squared(TLE_T3.4.lme_KL)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.06 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##SS
TLE_T3.4.lme_SS<-lmer(log(TLE_cm.day+1)~Origin+(1|Genotype), data=Growth_T3.4[which(Growth_T3.4$Site=="SS"),])
summary(TLE_T3.4.lme_SS)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(TLE_cm.day + 1) ~ Origin + (1 | Genotype)
Data: Growth_T3.4[which(Growth_T3.4$Site == "SS"), ]
REML criterion at convergence: -185.6
Scaled residuals:
Min 1Q Median 3Q Max
-2.44083 -0.56044 0.02132 0.52475 2.84423
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.001374 0.03707
Residual 0.001934 0.04398
Number of obs: 60, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.11023 0.02286 2.27122 4.822 0.03125 *
OriginTransplant 0.03439 0.01135 56.00000 3.029 0.00371 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.248
eta_squared(TLE_T3.4.lme_SS)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.14 | [0.03, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save results
TLE_T3.4.site.res<-data.frame(rbind(summary(TLE_T3.4.lme_KL)$coefficients[-1,],
summary(TLE_T3.4.lme_SS)$coefficients[-1,]))
TLE_T3.4.site.res$Predictor<-c("KL Origin", "SS Origin")
TLE_T3.4.site.res$EtaSq<-c(eta_squared(TLE_T3.4.lme_KL)$Eta2, eta_squared(TLE_T3.4.lme_SS)$Eta2)
TLE_T3.4.site.res$Response<-rep("Growth", nrow(TLE_T3.4.site.res))
TLE_T3.4.site.res$Timepoint<-rep("TP3_4", nrow(TLE_T3.4.site.res))
##Combine results
TLE_T3.4.res<-rbind(TLE_T3.4.res, TLE_T3.4.site.res)
##Summary statistics by Site and Origin
TP3.4_Growth.sum<-summarySE(Growth_T3.4, measurevar="TLE_cm.day", groupvars=c("Site", "Origin", "Site.Orig"), na.rm=TRUE)
##Plot Average Growth Rate across Treatments
TP3.4_Growth.plot<-ggplot(TP3.4_Growth.sum, aes(x=Site, y=TLE_cm.day, colour=Site.Orig)) +
scale_colour_manual(values=Orig.colors.o)+
geom_errorbar(aes(ymin=TLE_cm.day-se, ymax=TLE_cm.day+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Origin", y=expression(paste('Growth Rate (cm day'^-1*")")), colour="Origin")+
ylim(0, 0.25)+
annotate("text", x=1, y=0.24, label="-", size=sig.sz, fontface="bold")+
annotate("text", x=2, y=0.2, label="**", size=sig.sz, fontface="bold")+
geom_bracket(xmin=1, xmax=2, y.position=0.01, label.size=levels.sz, label="***", inherit.aes = FALSE); TP3.4_Growth.plot
SS
##Summary statistics by Genotype and Origin
TP3.4_Growth.sum.SS<-summarySE(subset(Growth_T3.4, Site=="SS"), measurevar="TLE_cm.day", groupvars=c("Genotype", "Origin"), na.rm=TRUE)
TP3.4_Growth.sum.SS
##Calculate Percent Difference between Native and Transplant by Genotype
TP3.4_Growth.SS<-TP3.4_Growth.sum.SS[which(TP3.4_Growth.sum.SS$Origin=="Native"), c(1,4)]
names(TP3.4_Growth.SS)[2]<-paste(names(TP3.4_Growth.SS)[2], "N", sep="_")
TP3.4_Growth.SS$TLE_cm.day_T<-TP3.4_Growth.sum.SS$TLE_cm.day[which(TP3.4_Growth.sum.SS$Origin=="Transplant")]
TP3.4_Growth.SS$Perc.Inc<-((TP3.4_Growth.SS$TLE_cm.day_T-TP3.4_Growth.SS$TLE_cm.day_N)/TP3.4_Growth.SS$TLE_cm.day_N)*100
TP3.4_Growth.SS
##Mean and Standard Error
mean(TP3.4_Growth.SS$Perc.Inc)
[1] 34.22395
std.error(TP3.4_Growth.SS$Perc.Inc)
[1] 18.61694
KL
##Summary statistics by Genotype and Origin
TP3.4_Growth.sum.KL<-summarySE(subset(Growth_T3.4, Site=="KL"), measurevar="TLE_cm.day", groupvars=c("Genotype", "Origin"), na.rm=TRUE)
TP3.4_Growth.sum.KL
##Calculate Percent Difference between Native and Transplant by Genotype
TP3.4_Growth.KL<-TP3.4_Growth.sum.KL[which(TP3.4_Growth.sum.KL$Origin=="Native"), c(1,4)]
names(TP3.4_Growth.KL)[2]<-paste(names(TP3.4_Growth.KL)[2], "N", sep="_")
TP3.4_Growth.KL$TLE_cm.day_T<-TP3.4_Growth.sum.KL$TLE_cm.day[which(TP3.4_Growth.sum.SS$Origin=="Transplant")]
TP3.4_Growth.KL$Perc.Inc<-((TP3.4_Growth.KL$TLE_cm.day_N-TP3.4_Growth.KL$TLE_cm.day_T)/TP3.4_Growth.KL$TLE_cm.day_T)*100
TP3.4_Growth.KL
##Mean and Standard Error
mean(TP3.4_Growth.KL$Perc.Inc)
[1] 21.37099
std.error(TP3.4_Growth.KL$Perc.Inc)
[1] 4.101183
##Growth Data dataframe
GrowthData<-Growth_T1.2
##Specify Timepoints
GrowthData<-GrowthData %>% dplyr::rename(T1.2_Ext_cm = Ext_cm, T1.2_TLE_cm.day = TLE_cm.day)
##Merge with Timepoints 3 to 4
GrowthData<-merge(GrowthData, Growth_T3.4, all=TRUE)
##Specify Timepoints
GrowthData<-GrowthData %>% dplyr::rename(T3.4_Ext_cm = Ext_cm, T3.4_TLE_cm.day = TLE_cm.day)
##Growth Data
write.csv(GrowthData, "Outputs/GrowthData.csv", row.names=FALSE)
##Control
Therm_C<-subset(Thermal, Treat=="C")
##Heated
Therm_H<-subset(Thermal, Treat=="H")
Calculating retention as proportion remaining relative to corresponding control levels (0-1) for each site and genotype at each timepoint.
##Calculate averages for Control Treatment for each Site, Genotype, and Timepoint
Therm_C<-na.omit(Therm_C)
names(Therm_C)
[1] "ID" "RandN" "TimeP" "Site" "Genotype" "Treat"
[7] "Treatment" "Orig" "Origin" "Set" "Site.Orig" "SA_cm2"
[13] "Chl_ug.cm2" "Sym10.6_cm2" "Fv_Fm"
Therm_C.a<-aggregate(Therm_C[,c(13:15)], list(Therm_C$Site, Therm_C$Genotype, Therm_C$TimeP, Therm_C$Origin), mean, na.action = na.omit)
names(Therm_C.a)[1:4]<-c("Site", "Genotype", "TimeP", "Origin")
names(Therm_C.a)[5:7]<-paste(names(Therm_C.a)[5:7], "C", sep="_")
##Merge Control Averages with Heated Samples
#Merges by Site, Genotype, and Timepoint
Therm_H<-merge(Therm_H, Therm_C.a, all.x=TRUE )
##Calculate Proportion Retained for each Bleaching Metric relative to Control
Therm_H$Chl.prop<-round((Therm_H$Chl_ug.cm2/Therm_H$Chl_ug.cm2_C), 4)
Therm_H$Sym.prop<-round((Therm_H$Sym10.6_cm2/Therm_H$Sym10.6_cm2_C), 4)
Therm_H$PAM.prop<-round((Therm_H$Fv_Fm/Therm_H$Fv_Fm_C), 4)
##Set values >1 to 1
Therm_H$Chl.prop[which(Therm_H$Chl.prop>1)]<-1.0000
Therm_H$Sym.prop[which(Therm_H$Sym.prop>1)]<-1.0000
Therm_H$PAM.prop[which(Therm_H$PAM.prop>1)]<-1.0000
##Subset Initial Timepoint
Thermal_IN<-subset(Thermal, TimeP=="IN")
Thermal_IN$Site.Treat<-paste(Thermal_IN$Site, Thermal_IN$Treat, sep=".")
Thermal_IN$Site.Treat<-factor(Thermal_IN$Site.Treat, levels=c("KL.C", "KL.H", "SS.C", "SS.H"))
Therm_IN<-subset(Therm_H, TimeP=="IN")
#Check normality
hist(Thermal_IN$Sym10.6_cm2)
shapiro.test(Thermal_IN$Sym10.6_cm2)
Shapiro-Wilk normality test
data: Thermal_IN$Sym10.6_cm2
W = 0.90824, p-value = 0.001334
#Not normal
hist(log(Thermal_IN$Sym10.6_cm2+1))
shapiro.test(log(Thermal_IN$Sym10.6_cm2+1))
Shapiro-Wilk normality test
data: log(Thermal_IN$Sym10.6_cm2 + 1)
W = 0.95667, p-value = 0.07936
#Normal
##Model
#Function of Site and Treatment, with Genotype as a Random effect
Sym.lme_IN<-lmer(log(Sym10.6_cm2+1)~Site*Treatment+(1|Genotype), data=Thermal_IN)
##Check residuals
Sym.lme_res_IN <- simulateResiduals(fittedModel = Sym.lme_IN, plot = F)
plot(Sym.lme_res_IN)
##Model results
summary(Sym.lme_IN)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Sym10.6_cm2 + 1) ~ Site * Treatment + (1 | Genotype)
Data: Thermal_IN
REML criterion at convergence: -29.3
Scaled residuals:
Min 1Q Median 3Q Max
-1.79325 -0.66013 0.08255 0.57806 2.12486
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.0008222 0.02867
Residual 0.0230525 0.15183
Number of obs: 47, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.38823 0.04685 14.06727 8.286 8.77e-07 ***
SiteSS 0.19508 0.06198 41.00308 3.147 0.00307 **
TreatmentHeat -0.20370 0.06340 41.10985 -3.213 0.00256 **
SiteSS:TreatmentHeat -0.01724 0.08867 41.05799 -0.194 0.84677
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) SiteSS TrtmnH
SiteSS -0.661
TreatmentHt -0.647 0.489
StSS:TrtmnH 0.462 -0.699 -0.715
eta_squared(Sym.lme_IN)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
----------------------------------------------
Site | 0.30 | [0.12, 1.00]
Treatment | 0.36 | [0.17, 1.00]
Site:Treatment | 9.20e-04 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
emmeans(Sym.lme_IN, pairwise ~ Site | Treatment)
$emmeans
Treatment = Control:
Site emmean SE df lower.CL upper.CL
KL 0.388 0.0469 14.1 0.2878 0.489
SS 0.583 0.0469 14.1 0.4829 0.684
Treatment = Heat:
Site emmean SE df lower.CL upper.CL
KL 0.185 0.0488 15.6 0.0808 0.288
SS 0.362 0.0469 14.1 0.2619 0.463
Degrees-of-freedom method: kenward-roger
Results are given on the log(mu + 1) (not the response) scale.
Confidence level used: 0.95
$contrasts
Treatment = Control:
contrast estimate SE df t.ratio p.value
KL - SS -0.195 0.0620 41.0 -3.147 0.0031
Treatment = Heat:
contrast estimate SE df t.ratio p.value
KL - SS -0.178 0.0635 41.1 -2.801 0.0077
Note: contrasts are still on the log(mu + 1) scale
Degrees-of-freedom method: kenward-roger
emmeans(Sym.lme_IN, pairwise ~ Treatment | Site)
$emmeans
Site = KL:
Treatment emmean SE df lower.CL upper.CL
Control 0.388 0.0469 14.1 0.2878 0.489
Heat 0.185 0.0488 15.6 0.0808 0.288
Site = SS:
Treatment emmean SE df lower.CL upper.CL
Control 0.583 0.0469 14.1 0.4829 0.684
Heat 0.362 0.0469 14.1 0.2619 0.463
Degrees-of-freedom method: kenward-roger
Results are given on the log(mu + 1) (not the response) scale.
Confidence level used: 0.95
$contrasts
Site = KL:
contrast estimate SE df t.ratio p.value
Control - Heat 0.204 0.0635 41.1 3.209 0.0026
Site = SS:
contrast estimate SE df t.ratio p.value
Control - Heat 0.221 0.0620 41.0 3.564 0.0009
Note: contrasts are still on the log(mu + 1) scale
Degrees-of-freedom method: kenward-roger
##Save model results
Sym_IN.res<-data.frame(summary(Sym.lme_IN)$coefficients[-1,])
##Add pairwise
names(Sym_IN.res)<-names(data.frame(emmeans(Sym.lme_IN, pairwise ~ Site | Treatment)$contrasts))[-c(1:2)]
Sym_IN.res<-rbind(Sym_IN.res,
data.frame(rbind(emmeans(Sym.lme_IN, pairwise ~ Site | Treatment)$contrasts[1], emmeans(Sym.lme_IN, pairwise ~ Site | Treatment)$contrasts[2], emmeans(Sym.lme_IN, pairwise ~ Treatment | Site)$contrasts[1], emmeans(Sym.lme_IN, pairwise ~ Treatment | Site)$contrasts[2]))[,-c(1:3)])
##Metadata
Sym_IN.res$Predictor<-c("Site", "Treatment", "Site x Treatment", "Control Site", "Heated Site", "KL Treatment", "SS Treatment")
Sym_IN.res$Response<-rep("Symbionts", nrow(Sym_IN.res))
Sym_IN.res$EtaSq<-c(eta_squared(Sym.lme_IN)$Eta2, rep(NA, 4))
##Summary statistics by Site
IN_Sym.sum<-summarySE(Thermal_IN, measurevar="Sym10.6_cm2", groupvars=c("Site", "Treatment", "Site.Treat"), na.rm=TRUE)
##Plot Average Symbionts across Treatments
IN_Sym.plot<-ggplot(IN_Sym.sum, aes(x=Site.Treat, y=Sym10.6_cm2, colour=Site)) +
scale_colour_manual(values=Site.colors.o)+
geom_errorbar(aes(ymin=Sym10.6_cm2-se, ymax=Sym10.6_cm2+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(aes(shape=Treatment), size=point.sz, position=position_dodge(width=0.5))+
scale_shape_manual(values=c(1,16))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Treatment", y=expression(paste('Symbiont Density ('*10^6,'cells cm'^-2*")")), colour="Site")+
ylim(0, 1); IN_Sym.plot
#Check normality
hist(Thermal_IN$Chl_ug.cm2)
shapiro.test(Thermal_IN$Chl_ug.cm2)
Shapiro-Wilk normality test
data: Thermal_IN$Chl_ug.cm2
W = 0.88059, p-value = 0.0002132
#Not normal
hist(log(Thermal_IN$Chl_ug.cm2+1))
shapiro.test(log(Thermal_IN$Chl_ug.cm2+1))
Shapiro-Wilk normality test
data: log(Thermal_IN$Chl_ug.cm2 + 1)
W = 0.92511, p-value = 0.005653
#Still not normal but less skewed
##Model
#Function of Site and Treatment, with Genotype as a Random effect
Chl.lme_IN<-lmer(log(Chl_ug.cm2+1)~Site*Treatment+(1|Genotype), data=Thermal_IN)
##Check residuals
Chl.lme_res_IN <- simulateResiduals(fittedModel = Chl.lme_IN, plot = F)
plot(Chl.lme_res_IN)
##Model results
summary(Chl.lme_IN)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Chl_ug.cm2 + 1) ~ Site * Treatment + (1 | Genotype)
Data: Thermal_IN
REML criterion at convergence: -69.1
Scaled residuals:
Min 1Q Median 3Q Max
-2.50589 -0.51517 0.07352 0.50169 2.35158
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.006137 0.07834
Residual 0.007930 0.08905
Number of obs: 46, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.45648 0.05202 2.95779 8.774 0.00329 **
SiteSS 0.27810 0.03635 39.99255 7.650 2.38e-09 ***
TreatmentHeat -0.37295 0.03721 40.00578 -10.024 1.80e-12 ***
SiteSS:TreatmentHeat -0.16117 0.05265 40.01278 -3.061 0.00393 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) SiteSS TrtmnH
SiteSS -0.349
TreatmentHt -0.341 0.489
StSS:TrtmnH 0.241 -0.691 -0.707
eta_squared(Chl.lme_IN)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
----------------------------------------------
Site | 0.58 | [0.41, 1.00]
Treatment | 0.88 | [0.82, 1.00]
Site:Treatment | 0.19 | [0.04, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
emmeans(Chl.lme_IN, pairwise ~ Site | Treatment)
$emmeans
Treatment = Control:
Site emmean SE df lower.CL upper.CL
KL 0.4565 0.0520 2.97 0.2899 0.623
SS 0.7346 0.0520 2.97 0.5680 0.901
Treatment = Heat:
Site emmean SE df lower.CL upper.CL
KL 0.0835 0.0526 3.11 -0.0808 0.248
SS 0.2005 0.0526 3.11 0.0361 0.365
Degrees-of-freedom method: kenward-roger
Results are given on the log(mu + 1) (not the response) scale.
Confidence level used: 0.95
$contrasts
Treatment = Control:
contrast estimate SE df t.ratio p.value
KL - SS -0.278 0.0364 40 -7.650 <.0001
Treatment = Heat:
contrast estimate SE df t.ratio p.value
KL - SS -0.117 0.0381 40 -3.070 0.0038
Note: contrasts are still on the log(mu + 1) scale
Degrees-of-freedom method: kenward-roger
emmeans(Chl.lme_IN, pairwise ~ Treatment | Site)
$emmeans
Site = KL:
Treatment emmean SE df lower.CL upper.CL
Control 0.4565 0.0520 2.97 0.2899 0.623
Heat 0.0835 0.0526 3.11 -0.0808 0.248
Site = SS:
Treatment emmean SE df lower.CL upper.CL
Control 0.7346 0.0520 2.97 0.5680 0.901
Heat 0.2005 0.0526 3.11 0.0361 0.365
Degrees-of-freedom method: kenward-roger
Results are given on the log(mu + 1) (not the response) scale.
Confidence level used: 0.95
$contrasts
Site = KL:
contrast estimate SE df t.ratio p.value
Control - Heat 0.373 0.0372 40 10.022 <.0001
Site = SS:
contrast estimate SE df t.ratio p.value
Control - Heat 0.534 0.0372 40 14.353 <.0001
Note: contrasts are still on the log(mu + 1) scale
Degrees-of-freedom method: kenward-roger
##Save model results
Chl_IN.res<-data.frame(summary(Chl.lme_IN)$coefficients[-1,])
##Add pairwise
names(Chl_IN.res)<-names(data.frame(emmeans(Chl.lme_IN, pairwise ~ Site | Treatment)$contrasts))[-c(1:2)]
Chl_IN.res<-rbind(Chl_IN.res,
data.frame(rbind(emmeans(Chl.lme_IN, pairwise ~ Site | Treatment)$contrasts[1], emmeans(Chl.lme_IN, pairwise ~ Site | Treatment)$contrasts[2], emmeans(Chl.lme_IN, pairwise ~ Treatment | Site)$contrasts[1], emmeans(Chl.lme_IN, pairwise ~ Treatment | Site)$contrasts[2]))[,-c(1:3)])
##Metadata
Chl_IN.res$Predictor<-c("Site", "Treatment", "Site x Treatment", "Control Site", "Heated Site", "KL Treatment", "SS Treatment")
Chl_IN.res$Response<-rep("Chlorophyll", nrow(Chl_IN.res))
Chl_IN.res$EtaSq<-c(eta_squared(Chl.lme_IN)$Eta2, rep(NA, 4))
##Summary statistics by Site
IN_Chl.sum<-summarySE(Thermal_IN, measurevar="Chl_ug.cm2", groupvars=c("Site", "Treatment", "Site.Treat"), na.rm=TRUE)
##Plot Average Chlorophyll across Treatments
IN_Chl.plot<-ggplot(IN_Chl.sum, aes(x=Site.Treat, y=Chl_ug.cm2, colour=Site)) +
scale_colour_manual(values=Site.colors.o)+
geom_errorbar(aes(ymin=Chl_ug.cm2-se, ymax=Chl_ug.cm2+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(aes(shape=Treatment), size=point.sz, position=position_dodge(width=0.5))+
scale_shape_manual(values=c(1,16))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Treatment", y=expression(paste('Total Chlorophyll (\u03BCg cm'^-2*")")), colour="Site")+
ylim(0, 1.25); IN_Chl.plot
#Check normality
hist(Thermal_IN$Fv_Fm)
shapiro.test(Thermal_IN$Fv_Fm)
Shapiro-Wilk normality test
data: Thermal_IN$Fv_Fm
W = 0.87006, p-value = 9.097e-05
#Not normal
hist(log(Thermal_IN$Fv_Fm+1))
shapiro.test(log(Thermal_IN$Fv_Fm+1))
Shapiro-Wilk normality test
data: log(Thermal_IN$Fv_Fm + 1)
W = 0.85365, p-value = 3.226e-05
#Still not normal but less skewed
##Model
#Function of Site and Treatment, with Genotype as a Random effect
PAM.lme_IN<-lmer(log(Fv_Fm+1)~Site*Treatment+(1|Genotype), data=Thermal_IN)
boundary (singular) fit: see help('isSingular')
##Check residuals
PAM.lme_res_IN <- simulateResiduals(fittedModel = PAM.lme_IN, plot = F)
plot(PAM.lme_res_IN)
##Model results
summary(PAM.lme_IN)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Fv_Fm + 1) ~ Site * Treatment + (1 | Genotype)
Data: Thermal_IN
REML criterion at convergence: -202.7
Scaled residuals:
Min 1Q Median 3Q Max
-3.7092 -0.4159 0.1391 0.6366 1.5720
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.0000000 0.00000
Residual 0.0004174 0.02043
Number of obs: 47, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.480202 0.005897 43.000000 81.425 < 2e-16 ***
SiteSS 0.013981 0.008340 43.000000 1.676 0.101
TreatmentHeat -0.044190 0.008528 43.000000 -5.182 5.57e-06 ***
SiteSS:TreatmentHeat 0.022162 0.011928 43.000000 1.858 0.070 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) SiteSS TrtmnH
SiteSS -0.707
TreatmentHt -0.692 0.489
StSS:TrtmnH 0.494 -0.699 -0.715
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')
eta_squared(PAM.lme_IN)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
----------------------------------------------
Site | 0.29 | [0.12, 1.00]
Treatment | 0.42 | [0.23, 1.00]
Site:Treatment | 0.07 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
emmeans(PAM.lme_IN, pairwise ~ Site | Treatment)
$emmeans
Treatment = Control:
Site emmean SE df lower.CL upper.CL
KL 0.480 0.00590 21.6 0.468 0.492
SS 0.494 0.00590 21.6 0.482 0.506
Treatment = Heat:
Site emmean SE df lower.CL upper.CL
KL 0.436 0.00618 23.5 0.423 0.449
SS 0.472 0.00590 21.6 0.460 0.484
Degrees-of-freedom method: kenward-roger
Results are given on the log(mu + 1) (not the response) scale.
Confidence level used: 0.95
$contrasts
Treatment = Control:
contrast estimate SE df t.ratio p.value
KL - SS -0.0140 0.00834 41.0 -1.676 0.1013
Treatment = Heat:
contrast estimate SE df t.ratio p.value
KL - SS -0.0361 0.00855 41.2 -4.229 0.0001
Note: contrasts are still on the log(mu + 1) scale
Degrees-of-freedom method: kenward-roger
emmeans(PAM.lme_IN, pairwise ~ Treatment | Site)
$emmeans
Site = KL:
Treatment emmean SE df lower.CL upper.CL
Control 0.480 0.00590 21.6 0.468 0.492
Heat 0.436 0.00618 23.5 0.423 0.449
Site = SS:
Treatment emmean SE df lower.CL upper.CL
Control 0.494 0.00590 21.6 0.482 0.506
Heat 0.472 0.00590 21.6 0.460 0.484
Degrees-of-freedom method: kenward-roger
Results are given on the log(mu + 1) (not the response) scale.
Confidence level used: 0.95
$contrasts
Site = KL:
contrast estimate SE df t.ratio p.value
Control - Heat 0.0442 0.00855 41.2 5.171 <.0001
Site = SS:
contrast estimate SE df t.ratio p.value
Control - Heat 0.0220 0.00834 41.0 2.641 0.0116
Note: contrasts are still on the log(mu + 1) scale
Degrees-of-freedom method: kenward-roger
##Save model results
PAM_IN.res<-data.frame(summary(PAM.lme_IN)$coefficients[-1,])
##Add pairwise
names(PAM_IN.res)<-names(data.frame(emmeans(PAM.lme_IN, pairwise ~ Site | Treatment)$contrasts))[-c(1:2)]
PAM_IN.res<-rbind(PAM_IN.res,
data.frame(rbind(emmeans(PAM.lme_IN, pairwise ~ Site | Treatment)$contrasts[1], emmeans(PAM.lme_IN, pairwise ~ Site | Treatment)$contrasts[2], emmeans(PAM.lme_IN, pairwise ~ Treatment | Site)$contrasts[1], emmeans(PAM.lme_IN, pairwise ~ Treatment | Site)$contrasts[2]))[,-c(1:3)])
##Metadata
PAM_IN.res$Predictor<-c("Site", "Treatment", "Site x Treatment", "Control Site", "Heated Site", "KL Treatment", "SS Treatment")
PAM_IN.res$Response<-rep("FvFm", nrow(PAM_IN.res))
PAM_IN.res$EtaSq<-c(eta_squared(PAM.lme_IN)$Eta2, rep(NA, 4))
Not normal, but model residuals are OK.
##Summary statistics by Site
IN_PAM.sum<-summarySE(Thermal_IN, measurevar="Fv_Fm", groupvars=c("Site", "Treatment", "Site.Treat"), na.rm=TRUE)
##Plot Average FvFm across Treatments
IN_PAM.plot<-ggplot(IN_PAM.sum, aes(x=Site.Treat, y=Fv_Fm, colour=Site)) +
scale_colour_manual(values=Site.colors.o)+
geom_errorbar(aes(ymin=Fv_Fm-se, ymax=Fv_Fm+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(aes(shape=Treatment), size=point.sz, position=position_dodge(width=0.5))+
scale_shape_manual(values=c(1,16))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Treatment", y="Photochemical Efficiency (Fv/Fm)", colour="Site")+
ylim(0.5, 0.655); IN_PAM.plot
#Check normality
hist(Therm_IN$Sym.prop)
shapiro.test(Therm_IN$Sym.prop)
Shapiro-Wilk normality test
data: Therm_IN$Sym.prop
W = 0.93042, p-value = 0.1117
#Normal
##Model
#Function of Site, with Genotype as a Random effect
Tol_Sym.lme_IN<-lmer(Sym.prop~Site+(1|Genotype), data=Therm_IN)
##Check residuals
Tol_Sym.lme_res_IN <- simulateResiduals(fittedModel = Tol_Sym.lme_IN, plot = F)
plot(Tol_Sym.lme_res_IN)
##Model results
summary(Tol_Sym.lme_IN)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: Sym.prop ~ Site + (1 | Genotype)
Data: Therm_IN
REML criterion at convergence: -2.8
Scaled residuals:
Min 1Q Median 3Q Max
-1.6802 -0.6578 -0.1521 0.6805 1.9852
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.07818 0.2796
Residual 0.03039 0.1743
Number of obs: 23, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.46139 0.16983 2.21175 2.717 0.1014
SiteSS 0.14316 0.07291 19.00328 1.964 0.0644 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
SiteSS -0.225
eta_squared(Tol_Sym.lme_IN)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Site | 0.17 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
Tol_Sym_IN.res<-data.frame(summary(Tol_Sym.lme_IN)$coefficients)
names(Tol_Sym_IN.res)<-names(Sym_IN.res)[1:5]
Tol_Sym_IN.res$Predictor<-rep("Site", nrow(Tol_Sym_IN.res))
Tol_Sym_IN.res$Response<-rep("Symbiont Retention", nrow(Tol_Sym_IN.res))
Tol_Sym_IN.res$EtaSq<-c(eta_squared(Tol_Sym.lme_IN)$Eta2)
##Combine results
Sym_IN.res<-rbind(Sym_IN.res, Tol_Sym_IN.res[2,])
##Summary statistics by Site
IN_TolSym.sum<-summarySE(Therm_IN, measurevar="Sym.prop", groupvars=c("Site"), na.rm=TRUE)
##Plot Average Symbiont Retention across Treatments
IN_TolSym.plot<-ggplot(IN_TolSym.sum, aes(x=Site, y=Sym.prop, colour=Site)) +
scale_colour_manual(values=Site.colors.o)+
geom_errorbar(aes(ymin=Sym.prop-se, ymax=Sym.prop+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site", y="Symbiont Retention", colour="Site")+
ylim(0, 1)+
annotate("text", x=1.5, y=0.75, label="-", size=sig.sz, fontface="bold"); IN_TolSym.plot
#Check normality
hist(Therm_IN$Chl.prop)
shapiro.test(Therm_IN$Chl.prop)
Shapiro-Wilk normality test
data: Therm_IN$Chl.prop
W = 0.98432, p-value = 0.9694
#Normal
##Model
#Function of Site, with Genotype as a Random effect
Tol_Chl.lme_IN<-lmer(Chl.prop~Site+(1|Genotype), data=Therm_IN)
##Check residuals
Tol_Chl.lme_res_IN <- simulateResiduals(fittedModel = Tol_Chl.lme_IN, plot = F)
plot(Tol_Chl.lme_res_IN)
##Model results
summary(Tol_Chl.lme_IN)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: Chl.prop ~ Site + (1 | Genotype)
Data: Therm_IN
REML criterion at convergence: -47.2
Scaled residuals:
Min 1Q Median 3Q Max
-1.3909 -0.5550 -0.1823 0.5284 1.9947
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.004576 0.06765
Residual 0.003436 0.05861
Number of obs: 22, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.16552 0.04290 2.40867 3.858 0.045 *
SiteSS 0.04306 0.02514 18.06310 1.713 0.104
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
SiteSS -0.293
eta_squared(Tol_Chl.lme_IN)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Site | 0.14 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
Tol_Chl_IN.res<-data.frame(summary(Tol_Chl.lme_IN)$coefficients)
names(Tol_Chl_IN.res)<-names(Chl_IN.res)[1:5]
Tol_Chl_IN.res$Predictor<-rep("Site", nrow(Tol_Chl_IN.res))
Tol_Chl_IN.res$Response<-rep("Chlorophyll Retention", nrow(Tol_Chl_IN.res))
Tol_Chl_IN.res$EtaSq<-c(eta_squared(Tol_Chl.lme_IN)$Eta2)
##Combine results
Chl_IN.res<-rbind(Chl_IN.res, Tol_Chl_IN.res[2,])
##Summary statistics by Site and Origin
IN_TolChl.sum<-summarySE(Therm_IN, measurevar="Chl.prop", groupvars=c("Site"), na.rm=TRUE)
##Plot Average Chlorophyll Retention across Treatments
IN_TolChl.plot<-ggplot(IN_TolChl.sum, aes(x=Site, y=Chl.prop, colour=Site)) +
scale_colour_manual(values=Site.colors.o)+
geom_errorbar(aes(ymin=Chl.prop-se, ymax=Chl.prop+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site", y="Chlorophyll Retention", colour="Site")+
ylim(0, 0.35); IN_TolChl.plot
#Check normality
hist(Therm_IN$PAM.prop)
shapiro.test(Therm_IN$PAM.prop)
Shapiro-Wilk normality test
data: Therm_IN$PAM.prop
W = 0.8289, p-value = 0.001154
#Not normal
hist(log(Therm_IN$PAM.prop+1))
shapiro.test(log(Therm_IN$PAM.prop+1))
Shapiro-Wilk normality test
data: log(Therm_IN$PAM.prop + 1)
W = 0.81061, p-value = 0.0005693
##Still not normal
##Model
#Function of Site, with Genotype as a Random effect
Tol_PAM.lme_IN<-lmer(PAM.prop~Site+(1|Genotype), data=Therm_IN)
##Check residuals
Tol_PAM.lme_res_IN <- simulateResiduals(fittedModel = Tol_PAM.lme_IN, plot = F)
plot(Tol_PAM.lme_res_IN)
##Model results
summary(Tol_PAM.lme_IN)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: PAM.prop ~ Site + (1 | Genotype)
Data: Therm_IN
REML criterion at convergence: -50.5
Scaled residuals:
Min 1Q Median 3Q Max
-2.7860 -0.2498 0.1225 0.5362 1.4591
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.002221 0.04713
Residual 0.003538 0.05948
Number of obs: 23, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.89138 0.03262 2.77551 27.325 0.000185 ***
SiteSS 0.05230 0.02487 19.00103 2.103 0.049053 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
SiteSS -0.399
eta_squared(Tol_PAM.lme_IN)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Site | 0.19 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
Tol_PAM_IN.res<-data.frame(summary(Tol_PAM.lme_IN)$coefficients)
names(Tol_PAM_IN.res)<-names(PAM_IN.res)[1:5]
Tol_PAM_IN.res$Predictor<-rep("Site", nrow(Tol_PAM_IN.res))
Tol_PAM_IN.res$Response<-rep("FvFm Retention", nrow(Tol_PAM_IN.res))
Tol_PAM_IN.res$EtaSq<-c(eta_squared(Tol_PAM.lme_IN)$Eta2)
##Combine results
PAM_IN.res<-rbind(PAM_IN.res, Tol_PAM_IN.res[2,])
Not normal, but model residuals are OK.
##Summary statistics by Site and Origin
IN_TolPAM.sum<-summarySE(Therm_IN, measurevar="PAM.prop", groupvars=c("Site"), na.rm=TRUE)
##Plot Average FvFm Retention across Treatments
IN_TolPAM.plot<-ggplot(IN_TolPAM.sum, aes(x=Site, y=PAM.prop, colour=Site)) +
scale_colour_manual(values=Site.colors.o)+
geom_errorbar(aes(ymin=PAM.prop-se, ymax=PAM.prop+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site", y="Photochemical Efficiency Retention", colour="Site")+
ylim(0.6, 1)+
annotate("text", x=1.5, y=1, label="*", size=sig.sz, fontface="bold"); IN_TolPAM.plot
##Subset Timepoint 1
Therm_TP1<-subset(Therm_H, TimeP=="TP1")
#Check normality
hist(Therm_TP1$Sym.prop)
shapiro.test(Therm_TP1$Sym.prop)
Shapiro-Wilk normality test
data: Therm_TP1$Sym.prop
W = 0.91699, p-value = 0.002332
#Not normal
hist(log(Therm_TP1$Sym.prop+1))
shapiro.test(log(Therm_TP1$Sym.prop+1))
Shapiro-Wilk normality test
data: log(Therm_TP1$Sym.prop + 1)
W = 0.92258, p-value = 0.003671
##Still not normal but less skewed
##Model with log +1 transformation
#Function of Site and Origin, with Genotype as a Random effect
#Interaction between Site and Origin
Tol_Sym.lme_TP1<-lmer(log(Sym.prop+1)~Origin*Site+(1|Genotype), data=Therm_TP1)
##Check residuals
Tol_Sym.lme_res_TP1 <- simulateResiduals(fittedModel = Tol_Sym.lme_TP1, plot = F)
plot(Tol_Sym.lme_res_TP1)
##Model results
summary(Tol_Sym.lme_TP1)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Sym.prop + 1) ~ Origin * Site + (1 | Genotype)
Data: Therm_TP1
REML criterion at convergence: -33.1
Scaled residuals:
Min 1Q Median 3Q Max
-3.1715 -0.5714 0.0059 0.7146 1.9799
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.03376 0.1837
Residual 0.01886 0.1373
Number of obs: 48, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.27621 0.11324 2.42412 2.439 0.1128
OriginTransplant 0.11668 0.05606 42.00000 2.081 0.0435 *
SiteSS 0.06309 0.05606 42.00000 1.125 0.2668
OriginTransplant:SiteSS -0.13250 0.07928 42.00000 -1.671 0.1021
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) OrgnTr SiteSS
OrgnTrnspln -0.248
SiteSS -0.248 0.500
OrgnTrn:SSS 0.175 -0.707 -0.707
eta_squared(Tol_Sym.lme_TP1)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-------------------------------------------
Origin | 0.04 | [0.00, 1.00]
Site | 1.51e-04 | [0.00, 1.00]
Origin:Site | 0.06 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
Tol_Sym_TP1.res<-data.frame(summary(Tol_Sym.lme_TP1)$coefficients[-1,])
Tol_Sym_TP1.res$Predictor<-c("Origin", "Site", "Origin x Site")
Tol_Sym_TP1.res$EtaSq<-c(eta_squared(Tol_Sym.lme_TP1)$Eta2)
Not normal, but model residuals are OK.
Effect size of Origin for each Site
##KL
Tol_Sym.lme_TP1_KL<-lmer(log(Sym.prop+1)~Origin+(1|Genotype), data=Therm_TP1[which(Therm_TP1$Site=="KL"),])
summary(Tol_Sym.lme_TP1_KL)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Sym.prop + 1) ~ Origin + (1 | Genotype)
Data: Therm_TP1[which(Therm_TP1$Site == "KL"), ]
REML criterion at convergence: -21.7
Scaled residuals:
Min 1Q Median 3Q Max
-2.22211 -0.64466 0.06469 0.75749 1.28798
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.04259 0.2064
Residual 0.01291 0.1136
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.27621 0.12358 2.14836 2.235 0.1461
OriginTransplant 0.11668 0.04638 20.00000 2.516 0.0205 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.188
eta_squared(Tol_Sym.lme_TP1_KL)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.24 | [0.02, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##SS
Tol_Sym.lme_TP1_SS<-lmer(log(Sym.prop+1)~Origin+(1|Genotype), data=Therm_TP1[which(Therm_TP1$Site=="SS"),])
summary(Tol_Sym.lme_TP1_SS)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Sym.prop + 1) ~ Origin + (1 | Genotype)
Data: Therm_TP1[which(Therm_TP1$Site == "SS"), ]
REML criterion at convergence: -8.8
Scaled residuals:
Min 1Q Median 3Q Max
-2.4670 -0.5662 0.1529 0.6190 1.4817
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.02355 0.1535
Residual 0.02580 0.1606
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.33930 0.10000 2.50724 3.393 0.0557 .
OriginTransplant -0.01583 0.06558 20.00000 -0.241 0.8118
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.328
eta_squared(Tol_Sym.lme_TP1_SS)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 2.90e-03 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save results
Tol_Sym_TP1.site.res<-data.frame(rbind(summary(Tol_Sym.lme_TP1_KL)$coefficients[-1,],
summary(Tol_Sym.lme_TP1_SS)$coefficients[-1,]))
Tol_Sym_TP1.site.res$Predictor<-c("KL Origin", "SS Origin")
Tol_Sym_TP1.site.res$EtaSq<-c(eta_squared(Tol_Sym.lme_TP1_KL)$Eta2, eta_squared(Tol_Sym.lme_TP1_SS)$Eta2)
##Combine results
Tol_Sym_TP1.res<-rbind(Tol_Sym_TP1.res, Tol_Sym_TP1.site.res)
Tol_Sym_TP1.res$Response<-rep("Symbiont Retention", nrow(Tol_Sym_TP1.res))
##Summary statistics by Site and Origin
TP1_TolSym.sum<-summarySE(Therm_TP1, measurevar="Sym.prop", groupvars=c("Site", "Origin", "Site.Orig"), na.rm=TRUE)
##Plot Average Symbiont Retention across Treatments
TP1_TolSym.plot<-ggplot(TP1_TolSym.sum, aes(x=Site, y=Sym.prop, colour=Site.Orig)) +
scale_colour_manual(values=Orig.colors.o)+
geom_errorbar(aes(ymin=Sym.prop-se, ymax=Sym.prop+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Origin", y="Symbiont Retention", colour="Origin")+
ylim(0, 1)+
annotate("text", x=1, y=0.75, label="*", size=sig.sz, fontface="bold"); TP1_TolSym.plot
#Check normality
hist(Therm_TP1$Chl.prop)
shapiro.test(Therm_TP1$Chl.prop)
Shapiro-Wilk normality test
data: Therm_TP1$Chl.prop
W = 0.94139, p-value = 0.01836
#Not normal
hist(log(Therm_TP1$Chl.prop+1))
shapiro.test(log(Therm_TP1$Chl.prop+1))
Shapiro-Wilk normality test
data: log(Therm_TP1$Chl.prop + 1)
W = 0.94765, p-value = 0.03227
##Still not normal but less skewed
##Model with log +1 transformation
#Function of Site and Origin, with Genotype as a Random effect
#Interaction between Site and Origin
Tol_Chl.lme_TP1<-lmer(log(Chl.prop+1)~Origin*Site+(1|Genotype), data=Therm_TP1)
##Check residuals
Tol_Chl.lme_res_TP1 <- simulateResiduals(fittedModel = Tol_Chl.lme_TP1, plot = F)
plot(Tol_Chl.lme_res_TP1)
##Model results
summary(Tol_Chl.lme_TP1)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Chl.prop + 1) ~ Origin * Site + (1 | Genotype)
Data: Therm_TP1
REML criterion at convergence: -112.2
Scaled residuals:
Min 1Q Median 3Q Max
-2.1439 -0.5458 -0.0382 0.6267 1.8760
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.008382 0.09155
Residual 0.003065 0.05537
Number of obs: 48, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.13931 0.05522 2.27664 2.523 0.113
OriginTransplant 0.01521 0.02260 42.00000 0.673 0.505
SiteSS 0.02628 0.02260 42.00000 1.163 0.251
OriginTransplant:SiteSS -0.02689 0.03197 42.00000 -0.841 0.405
Correlation of Fixed Effects:
(Intr) OrgnTr SiteSS
OrgnTrnspln -0.205
SiteSS -0.205 0.500
OrgnTrn:SSS 0.145 -0.707 -0.707
eta_squared(Tol_Chl.lme_TP1)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-------------------------------------------
Origin | 2.90e-04 | [0.00, 1.00]
Site | 0.02 | [0.00, 1.00]
Origin:Site | 0.02 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
Tol_Chl_TP1.res<-data.frame(summary(Tol_Chl.lme_TP1)$coefficients[-1,])
Tol_Chl_TP1.res$Predictor<-c("Origin", "Site", "Origin x Site")
Tol_Chl_TP1.res$EtaSq<-c(eta_squared(Tol_Chl.lme_TP1)$Eta2)
Effect size of Origin for each Site
##KL
Tol_Chl.lme_TP1_KL<-lmer(log(Chl.prop+1)~Origin+(1|Genotype), data=Therm_TP1[which(Therm_TP1$Site=="KL"),])
summary(Tol_Chl.lme_TP1_KL)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Chl.prop + 1) ~ Origin + (1 | Genotype)
Data: Therm_TP1[which(Therm_TP1$Site == "KL"), ]
REML criterion at convergence: -55.4
Scaled residuals:
Min 1Q Median 3Q Max
-2.04878 -0.58286 0.01887 0.57011 1.44796
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.007318 0.08555
Residual 0.002846 0.05334
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.13931 0.05173 2.18924 2.693 0.104
OriginTransplant 0.01521 0.02178 20.00000 0.698 0.493
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.210
eta_squared(Tol_Chl.lme_TP1_KL)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.02 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##SS
Tol_Chl.lme_TP1_SS<-lmer(log(Chl.prop+1)~Origin+(1|Genotype), data=Therm_TP1[which(Therm_TP1$Site=="SS"),])
summary(Tol_Chl.lme_TP1_SS)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Chl.prop + 1) ~ Origin + (1 | Genotype)
Data: Therm_TP1[which(Therm_TP1$Site == "SS"), ]
REML criterion at convergence: -50.8
Scaled residuals:
Min 1Q Median 3Q Max
-1.43149 -0.56754 -0.06902 0.58853 1.85344
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.009135 0.09558
Residual 0.003511 0.05926
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.16560 0.05777 2.18712 2.866 0.0931 .
OriginTransplant -0.01168 0.02419 20.00000 -0.483 0.6344
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.209
eta_squared(Tol_Chl.lme_TP1_SS)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.01 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save results
Tol_Chl_TP1.site.res<-data.frame(rbind(summary(Tol_Chl.lme_TP1_KL)$coefficients[-1,],
summary(Tol_Chl.lme_TP1_SS)$coefficients[-1,]))
Tol_Chl_TP1.site.res$Predictor<-c("KL Origin", "SS Origin")
Tol_Chl_TP1.site.res$EtaSq<-c(eta_squared(Tol_Chl.lme_TP1_KL)$Eta2, eta_squared(Tol_Chl.lme_TP1_SS)$Eta2)
##Combine results
Tol_Chl_TP1.res<-rbind(Tol_Chl_TP1.res, Tol_Chl_TP1.site.res)
Tol_Chl_TP1.res$Response<-rep("Chlorophyll Retention", nrow(Tol_Chl_TP1.res))
##Summary statistics by Site and Origin
TP1_TolChl.sum<-summarySE(Therm_TP1, measurevar="Chl.prop", groupvars=c("Site", "Origin", "Site.Orig"), na.rm=TRUE)
##Plot Average Chlorophyll Retention across Treatments
TP1_TolChl.plot<-ggplot(TP1_TolChl.sum, aes(x=Site, y=Chl.prop, colour=Site.Orig)) +
scale_colour_manual(values=Orig.colors.o)+
geom_errorbar(aes(ymin=Chl.prop-se, ymax=Chl.prop+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Origin", y="Chlorophyll Retention", colour="Origin")+
ylim(0, 0.35); TP1_TolChl.plot
#Check normality
hist(Therm_TP1$PAM.prop)
shapiro.test(Therm_TP1$PAM.prop)
Shapiro-Wilk normality test
data: Therm_TP1$PAM.prop
W = 0.84073, p-value = 1.245e-05
#Not normal
hist(log(Therm_TP1$PAM.prop+1))
shapiro.test(log(Therm_TP1$PAM.prop+1))
Shapiro-Wilk normality test
data: log(Therm_TP1$PAM.prop + 1)
W = 0.82049, p-value = 3.908e-06
##Still not normal
##Model
#Function of Site and Origin, with Genotype as a Random effect
#Interaction between Site and Origin
Tol_PAM.lme_TP1<-lmer(PAM.prop~Origin*Site+(1|Genotype), data=Therm_TP1)
##Check residuals
Tol_PAM.lme_res_TP1 <- simulateResiduals(fittedModel = Tol_PAM.lme_TP1, plot = F)
plot(Tol_PAM.lme_res_TP1)
##Model results
summary(Tol_PAM.lme_TP1)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: PAM.prop ~ Origin * Site + (1 | Genotype)
Data: Therm_TP1
REML criterion at convergence: -84.6
Scaled residuals:
Min 1Q Median 3Q Max
-2.3058 -0.3498 0.1733 0.6717 1.7818
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.002735 0.05230
Residual 0.006214 0.07883
Number of obs: 48, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.86181 0.03781 3.74555 22.794 3.72e-05 ***
OriginTransplant 0.06813 0.03218 42.00000 2.117 0.0402 *
SiteSS -0.01277 0.03218 42.00000 -0.397 0.6936
OriginTransplant:SiteSS -0.07942 0.04551 42.00000 -1.745 0.0883 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) OrgnTr SiteSS
OrgnTrnspln -0.426
SiteSS -0.426 0.500
OrgnTrn:SSS 0.301 -0.707 -0.707
eta_squared(Tol_PAM.lme_TP1)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-------------------------------------------
Origin | 0.04 | [0.00, 1.00]
Site | 0.11 | [0.01, 1.00]
Origin:Site | 0.07 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
Tol_PAM_TP1.res<-data.frame(summary(Tol_PAM.lme_TP1)$coefficients[-1,])
Tol_PAM_TP1.res$Predictor<-c("Origin", "Site", "Origin x Site")
Tol_PAM_TP1.res$EtaSq<-c(eta_squared(Tol_PAM.lme_TP1)$Eta2)
Not normal, but model residuals are OK.
Effect size of Origin for each Site
##KL
Tol_PAM.lme_TP1_KL<-lmer(PAM.prop~Origin+(1|Genotype), data=Therm_TP1[which(Therm_TP1$Site=="KL"),])
summary(Tol_PAM.lme_TP1_KL)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: PAM.prop ~ Origin + (1 | Genotype)
Data: Therm_TP1[which(Therm_TP1$Site == "KL"), ]
REML criterion at convergence: -57.2
Scaled residuals:
Min 1Q Median 3Q Max
-2.21445 -0.36299 0.09321 0.61667 1.51097
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.0009092 0.03015
Residual 0.0031116 0.05578
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.86181 0.02371 3.34791 36.341 1.82e-05 ***
OriginTransplant 0.06813 0.02277 20.00000 2.992 0.00721 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.480
eta_squared(Tol_PAM.lme_TP1_KL)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.31 | [0.06, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##SS
Tol_PAM.lme_TP1_SS<-lmer(PAM.prop~Origin+(1|Genotype), data=Therm_TP1[which(Therm_TP1$Site=="SS"),])
summary(Tol_PAM.lme_TP1_SS)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: PAM.prop ~ Origin + (1 | Genotype)
Data: Therm_TP1[which(Therm_TP1$Site == "SS"), ]
REML criterion at convergence: -33.3
Scaled residuals:
Min 1Q Median 3Q Max
-1.69510 -0.36658 0.09367 0.51233 1.74473
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.005306 0.07284
Residual 0.008775 0.09367
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.84904 0.05000 2.73586 16.981 0.000746 ***
OriginTransplant -0.01128 0.03824 20.00000 -0.295 0.770994
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.382
eta_squared(Tol_PAM.lme_TP1_SS)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 4.33e-03 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save results
Tol_PAM_TP1.site.res<-data.frame(rbind(summary(Tol_PAM.lme_TP1_KL)$coefficients[-1,],
summary(Tol_PAM.lme_TP1_SS)$coefficients[-1,]))
Tol_PAM_TP1.site.res$Predictor<-c("KL Origin", "SS Origin")
Tol_PAM_TP1.site.res$EtaSq<-c(eta_squared(Tol_PAM.lme_TP1_KL)$Eta2, eta_squared(Tol_PAM.lme_TP1_SS)$Eta2)
##Combine results
Tol_PAM_TP1.res<-rbind(Tol_PAM_TP1.res, Tol_PAM_TP1.site.res)
Tol_PAM_TP1.res$Response<-rep("FvFm Retention", nrow(Tol_PAM_TP1.res))
##Summary statistics by Site and Origin
TP1_TolPAM.sum<-summarySE(Therm_TP1, measurevar="PAM.prop", groupvars=c("Site", "Origin", "Site.Orig"), na.rm=TRUE)
##Plot Average FvFm Retention across Treatments
TP1_TolPAM.plot<-ggplot(TP1_TolPAM.sum, aes(x=Site, y=PAM.prop, colour=Site.Orig)) +
scale_colour_manual(values=Orig.colors.o)+
geom_errorbar(aes(ymin=PAM.prop-se, ymax=PAM.prop+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Origin", y="Photochemical Efficiency Retention", colour="Origin")+
ylim(0.6, 1)+
annotate("text", x=1, y=1, label="**", size=sig.sz, fontface="bold"); TP1_TolPAM.plot
##Combine Results
Tol_T1.res<-rbind(Tol_PAM_TP1.res, Tol_Sym_TP1.res, Tol_Chl_TP1.res)
##Add Timepoint
Tol_T1.res$TimeP<-rep("TP1", nrow(Tol_T1.res))
##Subset Timepoint 2
Therm_TP2<-subset(Therm_H, TimeP=="TP2")
#Check normality
hist(Therm_TP2$Sym.prop)
shapiro.test(Therm_TP2$Sym.prop)
Shapiro-Wilk normality test
data: Therm_TP2$Sym.prop
W = 0.89527, p-value = 0.0004439
#Not normal
hist(log(Therm_TP2$Sym.prop+1))
shapiro.test(log(Therm_TP2$Sym.prop+1))
Shapiro-Wilk normality test
data: log(Therm_TP2$Sym.prop + 1)
W = 0.88978, p-value = 0.000299
##Still not normal
##Model
#Function of Site and Origin, with Genotype as a Random effect
#Interaction between Site and Origin
Tol_Sym.lme_TP2<-lmer(Sym.prop~Origin*Site+(1|Genotype), data=Therm_TP2)
##Check residuals
Tol_Sym.lme_res_TP2 <- simulateResiduals(fittedModel = Tol_Sym.lme_TP2, plot = F)
plot(Tol_Sym.lme_res_TP2)
##Model results
summary(Tol_Sym.lme_TP2)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: Sym.prop ~ Origin * Site + (1 | Genotype)
Data: Therm_TP2
REML criterion at convergence: -3.4
Scaled residuals:
Min 1Q Median 3Q Max
-2.20381 -0.54227 0.01202 0.47665 2.09092
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.04681 0.2164
Residual 0.03767 0.1941
Number of obs: 48, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.805833 0.136905 2.613359 5.886 0.0142 *
OriginTransplant -0.107900 0.079237 42.000000 -1.362 0.1805
SiteSS -0.178425 0.079237 42.000000 -2.252 0.0296 *
OriginTransplant:SiteSS 0.005183 0.112058 42.000000 0.046 0.9633
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) OrgnTr SiteSS
OrgnTrnspln -0.289
SiteSS -0.289 0.500
OrgnTrn:SSS 0.205 -0.707 -0.707
eta_squared(Tol_Sym.lme_TP2)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-------------------------------------------
Origin | 0.08 | [0.00, 1.00]
Site | 0.19 | [0.04, 1.00]
Origin:Site | 5.09e-05 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
Tol_Sym_TP2.res<-data.frame(summary(Tol_Sym.lme_TP2)$coefficients[-1,])
Tol_Sym_TP2.res$Predictor<-c("Origin", "Site", "Origin x Site")
Tol_Sym_TP2.res$EtaSq<-c(eta_squared(Tol_Sym.lme_TP2)$Eta2)
Not normal, but model residuals are OK.
Effect size of Origin for each Site
##KL
Tol_Sym.lme_TP2_KL<-lmer(Sym.prop~Origin+(1|Genotype), data=Therm_TP2[which(Therm_TP2$Site=="KL"),])
summary(Tol_Sym.lme_TP2_KL)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: Sym.prop ~ Origin + (1 | Genotype)
Data: Therm_TP2[which(Therm_TP2$Site == "KL"), ]
REML criterion at convergence: -1.7
Scaled residuals:
Min 1Q Median 3Q Max
-2.0500 -0.8276 0.2039 0.4804 2.0614
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.04096 0.2024
Residual 0.03487 0.1867
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.80583 0.12868 2.40096 6.262 0.0154 *
OriginTransplant -0.10790 0.07623 20.00000 -1.415 0.1723
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.296
eta_squared(Tol_Sym.lme_TP2_KL)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.09 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##SS
Tol_Sym.lme_TP2_SS<-lmer(Sym.prop~Origin+(1|Genotype), data=Therm_TP2[which(Therm_TP2$Site=="SS"),])
summary(Tol_Sym.lme_TP2_SS)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: Sym.prop ~ Origin + (1 | Genotype)
Data: Therm_TP2[which(Therm_TP2$Site == "SS"), ]
REML criterion at convergence: 3.2
Scaled residuals:
Min 1Q Median 3Q Max
-1.98373 -0.46709 -0.03333 0.50742 1.82777
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.04792 0.2189
Residual 0.04392 0.2096
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.62741 0.14011 2.42980 4.478 0.032 *
OriginTransplant -0.10272 0.08556 20.00000 -1.201 0.244
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.305
eta_squared(Tol_Sym.lme_TP2_SS)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.07 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save results
Tol_Sym_TP2.site.res<-data.frame(rbind(summary(Tol_Sym.lme_TP2_KL)$coefficients[-1,],
summary(Tol_Sym.lme_TP2_SS)$coefficients[-1,]))
Tol_Sym_TP2.site.res$Predictor<-c("KL Origin", "SS Origin")
Tol_Sym_TP2.site.res$EtaSq<-c(eta_squared(Tol_Sym.lme_TP2_KL)$Eta2, eta_squared(Tol_Sym.lme_TP2_SS)$Eta2)
##Combine results
Tol_Sym_TP2.res<-rbind(Tol_Sym_TP2.res, Tol_Sym_TP2.site.res)
Tol_Sym_TP2.res$Response<-rep("Symbiont Retention", nrow(Tol_Sym_TP2.res))
##Summary statistics by Site and Origin
TP2_TolSym.sum<-summarySE(Therm_TP2, measurevar="Sym.prop", groupvars=c("Site", "Origin", "Site.Orig"), na.rm=TRUE)
##Plot Average Retention across Treatments
TP2_TolSym.plot<-ggplot(TP2_TolSym.sum, aes(x=Site, y=Sym.prop, colour=Site.Orig)) +
scale_colour_manual(values=Orig.colors.o)+
geom_errorbar(aes(ymin=Sym.prop-se, ymax=Sym.prop+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Origin", y="Symbiont Retention", colour="Origin")+
ylim(0, 1); TP2_TolSym.plot
#Check normality
hist(Therm_TP2$Chl.prop)
shapiro.test(Therm_TP2$Chl.prop)
Shapiro-Wilk normality test
data: Therm_TP2$Chl.prop
W = 0.97677, p-value = 0.4523
#Normal
##Model
#Function of Site and Origin, with Genotype as a Random effect
#Interaction between Site and Origin
Tol_Chl.lme_TP2<-lmer(Chl.prop~Origin*Site+(1|Genotype), data=Therm_TP2)
##Check residuals
Tol_Chl.lme_res_TP2 <- simulateResiduals(fittedModel = Tol_Chl.lme_TP2, plot = F)
plot(Tol_Chl.lme_res_TP2)
##Model results
summary(Tol_Chl.lme_TP2)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: Chl.prop ~ Origin * Site + (1 | Genotype)
Data: Therm_TP2
REML criterion at convergence: -81
Scaled residuals:
Min 1Q Median 3Q Max
-2.0024 -0.7622 -0.1174 0.6824 2.9666
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.006146 0.07840
Residual 0.006538 0.08086
Number of obs: 48, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.31077 0.05093 2.81331 6.102 0.0106 *
OriginTransplant -0.03844 0.03301 42.00000 -1.165 0.2508
SiteSS -0.06514 0.03301 42.00000 -1.973 0.0550 .
OriginTransplant:SiteSS 0.06166 0.04668 42.00000 1.321 0.1937
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) OrgnTr SiteSS
OrgnTrnspln -0.324
SiteSS -0.324 0.500
OrgnTrn:SSS 0.229 -0.707 -0.707
eta_squared(Tol_Chl.lme_TP2)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-------------------------------------------
Origin | 2.53e-03 | [0.00, 1.00]
Site | 0.05 | [0.00, 1.00]
Origin:Site | 0.04 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
Tol_Chl_TP2.res<-data.frame(summary(Tol_Chl.lme_TP2)$coefficients[-1,])
Tol_Chl_TP2.res$Predictor<-c("Origin", "Site", "Origin x Site")
Tol_Chl_TP2.res$EtaSq<-c(eta_squared(Tol_Chl.lme_TP2)$Eta2)
Effect size of Origin for each Site
##KL
Tol_Chl.lme_TP2_KL<-lmer(Chl.prop~Origin+(1|Genotype), data=Therm_TP2[which(Therm_TP2$Site=="KL"),])
summary(Tol_Chl.lme_TP2_KL)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: Chl.prop ~ Origin + (1 | Genotype)
Data: Therm_TP2[which(Therm_TP2$Site == "KL"), ]
REML criterion at convergence: -39.9
Scaled residuals:
Min 1Q Median 3Q Max
-0.9725 -0.7459 -0.3411 0.6270 2.6035
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.006040 0.07772
Residual 0.006269 0.07918
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.31077 0.05036 2.48246 6.171 0.0146 *
OriginTransplant -0.03844 0.03232 20.00001 -1.189 0.2483
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.321
eta_squared(Tol_Chl.lme_TP2_KL)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.07 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##SS
Tol_Chl.lme_TP2_SS<-lmer(Chl.prop~Origin+(1|Genotype), data=Therm_TP2[which(Therm_TP2$Site=="SS"),])
summary(Tol_Chl.lme_TP2_SS)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: Chl.prop ~ Origin + (1 | Genotype)
Data: Therm_TP2[which(Therm_TP2$Site == "SS"), ]
REML criterion at convergence: -46
Scaled residuals:
Min 1Q Median 3Q Max
-2.2119 -0.6418 0.2701 0.7660 1.5135
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.009505 0.09749
Residual 0.004441 0.06664
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.24563 0.05948 2.22611 4.129 0.0448 *
OriginTransplant 0.02322 0.02720 20.00000 0.853 0.4035
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.229
eta_squared(Tol_Chl.lme_TP2_SS)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.04 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save results
Tol_Chl_TP2.site.res<-data.frame(rbind(summary(Tol_Chl.lme_TP2_KL)$coefficients[-1,],
summary(Tol_Chl.lme_TP2_SS)$coefficients[-1,]))
Tol_Chl_TP2.site.res$Predictor<-c("KL Origin", "SS Origin")
Tol_Chl_TP2.site.res$EtaSq<-c(eta_squared(Tol_Chl.lme_TP2_KL)$Eta2, eta_squared(Tol_Chl.lme_TP2_SS)$Eta2)
##Combine results
Tol_Chl_TP2.res<-rbind(Tol_Chl_TP2.res, Tol_Chl_TP2.site.res)
Tol_Chl_TP2.res$Response<-rep("Chlorophyll Retention", nrow(Tol_Chl_TP2.res))
##Summary statistics by Site and Origin
TP2_TolChl.sum<-summarySE(Therm_TP2, measurevar="Chl.prop", groupvars=c("Site", "Origin", "Site.Orig"), na.rm=TRUE)
##Plot Average Retention across Treatments
TP2_TolChl.plot<-ggplot(TP2_TolChl.sum, aes(x=Site, y=Chl.prop, colour=Site.Orig)) +
scale_colour_manual(values=Orig.colors.o)+
geom_errorbar(aes(ymin=Chl.prop-se, ymax=Chl.prop+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Origin", y="Chlorophyll Retention", colour="Origin")+
ylim(0, 0.35); TP2_TolChl.plot
#Check normality
hist(Therm_TP2$PAM.prop)
shapiro.test(Therm_TP2$PAM.prop)
Shapiro-Wilk normality test
data: Therm_TP2$PAM.prop
W = 0.73527, p-value = 5.896e-08
#Not normal
hist(log(Therm_TP2$PAM.prop+1))
shapiro.test(log(Therm_TP2$PAM.prop+1))
Shapiro-Wilk normality test
data: log(Therm_TP2$PAM.prop + 1)
W = 0.70205, p-value = 1.443e-08
##Still not normal
##Model
#Function of Site and Origin, with Genotype as a Random effect
#Interaction between Site and Origin
Tol_PAM.lme_TP2<-lmer(PAM.prop~Origin*Site+(1|Genotype), data=Therm_TP2)
##Check residuals
Tol_PAM.lme_res_TP2 <- simulateResiduals(fittedModel = Tol_PAM.lme_TP2, plot = F)
plot(Tol_PAM.lme_res_TP2)
##Model results
summary(Tol_PAM.lme_TP2)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: PAM.prop ~ Origin * Site + (1 | Genotype)
Data: Therm_TP2
REML criterion at convergence: -102.1
Scaled residuals:
Min 1Q Median 3Q Max
-3.8326 -0.3753 -0.0042 0.5350 1.7496
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.001675 0.04093
Residual 0.004192 0.06475
Number of obs: 48, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.85622 0.03013 3.92146 28.419 1.09e-05 ***
OriginTransplant 0.07526 0.02643 42.00000 2.847 0.00680 **
SiteSS 0.06487 0.02643 42.00000 2.454 0.01834 *
OriginTransplant:SiteSS -0.10206 0.03738 42.00000 -2.730 0.00921 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) OrgnTr SiteSS
OrgnTrnspln -0.439
SiteSS -0.439 0.500
OrgnTrn:SSS 0.310 -0.707 -0.707
eta_squared(Tol_PAM.lme_TP2)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-------------------------------------------
Origin | 0.04 | [0.00, 1.00]
Site | 0.01 | [0.00, 1.00]
Origin:Site | 0.15 | [0.02, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
Tol_PAM_TP2.res<-data.frame(summary(Tol_PAM.lme_TP2)$coefficients[-1,])
Tol_PAM_TP2.res$Predictor<-c("Origin", "Site", "Origin x Site")
Tol_PAM_TP2.res$EtaSq<-c(eta_squared(Tol_PAM.lme_TP2)$Eta2)
Not normal, but model residuals are OK.
Effect size of Origin for each Site
##KL
Tol_PAM.lme_TP2_KL<-lmer(PAM.prop~Origin+(1|Genotype), data=Therm_TP2[which(Therm_TP2$Site=="KL"),])
summary(Tol_PAM.lme_TP2_KL)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: PAM.prop ~ Origin + (1 | Genotype)
Data: Therm_TP2[which(Therm_TP2$Site == "KL"), ]
REML criterion at convergence: -40.2
Scaled residuals:
Min 1Q Median 3Q Max
-2.9208 -0.4020 0.2520 0.5255 1.0060
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.003195 0.05653
Residual 0.006508 0.08067
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.85622 0.04009 2.88245 21.356 0.00029 ***
OriginTransplant 0.07526 0.03293 19.99994 2.285 0.03336 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.411
eta_squared(Tol_PAM.lme_TP2_KL)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.21 | [0.01, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##SS
Tol_PAM.lme_TP2_SS<-lmer(PAM.prop~Origin+(1|Genotype), data=Therm_TP2[which(Therm_TP2$Site=="SS"),])
summary(Tol_PAM.lme_TP2_SS)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: PAM.prop ~ Origin + (1 | Genotype)
Data: Therm_TP2[which(Therm_TP2$Site == "SS"), ]
REML criterion at convergence: -71.4
Scaled residuals:
Min 1Q Median 3Q Max
-2.04910 -0.74292 -0.04708 0.53594 2.26126
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.0004937 0.02222
Residual 0.0016296 0.04037
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.92109 0.01733 3.31062 53.146 5.72e-06 ***
OriginTransplant -0.02680 0.01648 20.00000 -1.626 0.12
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.475
eta_squared(Tol_PAM.lme_TP2_SS)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.12 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save results
Tol_PAM_TP2.site.res<-data.frame(rbind(summary(Tol_PAM.lme_TP2_KL)$coefficients[-1,],
summary(Tol_PAM.lme_TP2_SS)$coefficients[-1,]))
Tol_PAM_TP2.site.res$Predictor<-c("KL Origin", "SS Origin")
Tol_PAM_TP2.site.res$EtaSq<-c(eta_squared(Tol_PAM.lme_TP2_KL)$Eta2, eta_squared(Tol_PAM.lme_TP2_SS)$Eta2)
##Combine results
Tol_PAM_TP2.res<-rbind(Tol_PAM_TP2.res, Tol_PAM_TP2.site.res)
Tol_PAM_TP2.res$Response<-rep("FvFm Retention", nrow(Tol_PAM_TP2.res))
##Summary statistics by Site and Origin
TP2_TolPAM.sum<-summarySE(Therm_TP2, measurevar="PAM.prop", groupvars=c("Site", "Origin", "Site.Orig"), na.rm=TRUE)
##Plot Average Retention across Treatments
TP2_TolPAM.plot<-ggplot(TP2_TolPAM.sum, aes(x=Site, y=PAM.prop, colour=Site.Orig)) +
scale_colour_manual(values=Orig.colors.o)+
geom_errorbar(aes(ymin=PAM.prop-se, ymax=PAM.prop+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Origin", y="Photochemical Efficiency Retention", colour="Origin")+
ylim(0.6, 1)+
annotate("text", x=1, y=1, label="*", size=sig.sz, fontface="bold"); TP2_TolPAM.plot
##Combine Results
Tol_T2.res<-rbind(Tol_PAM_TP2.res, Tol_Sym_TP2.res, Tol_Chl_TP2.res)
##Add Timepoint
Tol_T2.res$TimeP<-rep("TP2", nrow(Tol_T2.res))
##Subset Timepoint 3
Therm_TP3<-subset(Therm_H, TimeP=="TP3")
#Check normality
hist(Therm_TP3$Sym.prop)
shapiro.test(Therm_TP3$Sym.prop)
Shapiro-Wilk normality test
data: Therm_TP3$Sym.prop
W = 0.78653, p-value = 6.531e-07
#Not normal
hist(log(Therm_TP3$Sym.prop+1))
shapiro.test(log(Therm_TP3$Sym.prop+1))
Shapiro-Wilk normality test
data: log(Therm_TP3$Sym.prop + 1)
W = 0.78571, p-value = 6.267e-07
##Still not normal but less skewed
##Model with log +1 transformation
#Function of Site and Origin, with Genotype as a Random effect
#Interaction between Site and Origin
Tol_Sym.lme_TP3<-lmer(log(Sym.prop+1)~Origin*Site+(1|Genotype), data=Therm_TP3)
##Check residuals
Tol_Sym.lme_res_TP3 <- simulateResiduals(fittedModel = Tol_Sym.lme_TP3, plot = F)
plot(Tol_Sym.lme_res_TP3)
##Model results
summary(Tol_Sym.lme_TP3)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Sym.prop + 1) ~ Origin * Site + (1 | Genotype)
Data: Therm_TP3
REML criterion at convergence: -108.7
Scaled residuals:
Min 1Q Median 3Q Max
-2.1559 -0.5389 0.3778 0.5963 1.6051
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 8.41e-05 0.00917
Residual 3.90e-03 0.06245
Number of obs: 48, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.59780 0.01879 16.87446 31.815 <2e-16 ***
OriginTransplant 0.03978 0.02550 42.00000 1.560 0.1262
SiteSS 0.05892 0.02550 42.00000 2.311 0.0258 *
OriginTransplant:SiteSS -0.03130 0.03606 42.00000 -0.868 0.3903
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) OrgnTr SiteSS
OrgnTrnspln -0.678
SiteSS -0.678 0.500
OrgnTrn:SSS 0.480 -0.707 -0.707
eta_squared(Tol_Sym.lme_TP3)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-------------------------------------------
Origin | 0.04 | [0.00, 1.00]
Site | 0.12 | [0.01, 1.00]
Origin:Site | 0.02 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
Tol_Sym_TP3.res<-data.frame(summary(Tol_Sym.lme_TP3)$coefficients[-1,])
Tol_Sym_TP3.res$Predictor<-c("Origin", "Site", "Origin x Site")
Tol_Sym_TP3.res$EtaSq<-c(eta_squared(Tol_Sym.lme_TP3)$Eta2)
Not normal, but model residuals are OK.
Effect size of Origin for each Site
##KL
Tol_Sym.lme_TP3_KL<-lmer(log(Sym.prop+1)~Origin+(1|Genotype), data=Therm_TP3[which(Therm_TP3$Site=="KL"),])
summary(Tol_Sym.lme_TP3_KL)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Sym.prop + 1) ~ Origin + (1 | Genotype)
Data: Therm_TP3[which(Therm_TP3$Site == "KL"), ]
REML criterion at convergence: -48.9
Scaled residuals:
Min 1Q Median 3Q Max
-1.66523 -0.53461 0.01898 0.71797 1.50619
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.0003833 0.01958
Residual 0.0048325 0.06952
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.59780 0.02303 5.00842 25.955 1.56e-06 ***
OriginTransplant 0.03978 0.02838 20.00000 1.402 0.176
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.616
eta_squared(Tol_Sym.lme_TP3_KL)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.09 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##SS
Tol_Sym.lme_TP3_SS<-lmer(log(Sym.prop+1)~Origin+(1|Genotype), data=Therm_TP3[which(Therm_TP3$Site=="SS"),])
boundary (singular) fit: see help('isSingular')
summary(Tol_Sym.lme_TP3_SS)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Sym.prop + 1) ~ Origin + (1 | Genotype)
Data: Therm_TP3[which(Therm_TP3$Site == "SS"), ]
REML criterion at convergence: -61.8
Scaled residuals:
Min 1Q Median 3Q Max
-2.6314 -0.4640 0.5270 0.5670 0.6869
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.000000 0.00000
Residual 0.002812 0.05302
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.656722 0.015307 22.000000 42.904 <2e-16 ***
OriginTransplant 0.008483 0.021647 22.000000 0.392 0.699
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.707
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')
eta_squared(Tol_Sym.lme_TP3_SS)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 6.93e-03 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save results
Tol_Sym_TP3.site.res<-data.frame(rbind(summary(Tol_Sym.lme_TP3_KL)$coefficients[-1,],
summary(Tol_Sym.lme_TP3_SS)$coefficients[-1,]))
Tol_Sym_TP3.site.res$Predictor<-c("KL Origin", "SS Origin")
Tol_Sym_TP3.site.res$EtaSq<-c(eta_squared(Tol_Sym.lme_TP3_KL)$Eta2, eta_squared(Tol_Sym.lme_TP3_SS)$Eta2)
##Combine results
Tol_Sym_TP3.res<-rbind(Tol_Sym_TP3.res, Tol_Sym_TP3.site.res)
Tol_Sym_TP3.res$Response<-rep("Symbiont Retention", nrow(Tol_Sym_TP3.res))
##Summary statistics by Site and Origin
TP3_TolSym.sum<-summarySE(Therm_TP3, measurevar="Sym.prop", groupvars=c("Site", "Origin", "Site.Orig"), na.rm=TRUE)
##Plot Average Retention across Treatments
TP3_TolSym.plot<-ggplot(TP3_TolSym.sum, aes(x=Site, y=Sym.prop, colour=Site.Orig)) +
scale_colour_manual(values=Orig.colors.o)+
geom_errorbar(aes(ymin=Sym.prop-se, ymax=Sym.prop+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Origin", y="Symbiont Retention", colour="Origin")+
ylim(0, 1); TP3_TolSym.plot
#Check normality
hist(Therm_TP3$Chl.prop)
shapiro.test(Therm_TP3$Chl.prop)
Shapiro-Wilk normality test
data: Therm_TP3$Chl.prop
W = 0.89534, p-value = 0.0004461
#Not normal
hist(log(Therm_TP3$Chl.prop+1))
shapiro.test(log(Therm_TP3$Chl.prop+1))
Shapiro-Wilk normality test
data: log(Therm_TP3$Chl.prop + 1)
W = 0.9211, p-value = 0.003251
##Still not normal but less skewed
##Model with log +1 transformation
#Function of Site and Origin, with Genotype as a Random effect
#Interaction between Site and Origin
Tol_Chl.lme_TP3<-lmer(log(Chl.prop+1)~Origin*Site+(1|Genotype), data=Therm_TP3)
boundary (singular) fit: see help('isSingular')
##Check residuals
Tol_Chl.lme_res_TP3 <- simulateResiduals(fittedModel = Tol_Chl.lme_TP3, plot = F)
plot(Tol_Chl.lme_res_TP3)
##Model results
summary(Tol_Chl.lme_TP3)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Chl.prop + 1) ~ Origin * Site + (1 | Genotype)
Data: Therm_TP3
REML criterion at convergence: -47.3
Scaled residuals:
Min 1Q Median 3Q Max
-1.7095 -0.6080 -0.1405 0.4999 1.7650
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.00000 0.0000
Residual 0.01594 0.1262
Number of obs: 48, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.43851 0.03644 44.00000 12.032 1.65e-15 ***
OriginTransplant 0.03181 0.05154 44.00000 0.617 0.540
SiteSS 0.03852 0.05154 44.00000 0.747 0.459
OriginTransplant:SiteSS -0.02511 0.07289 44.00000 -0.344 0.732
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) OrgnTr SiteSS
OrgnTrnspln -0.707
SiteSS -0.707 0.500
OrgnTrn:SSS 0.500 -0.707 -0.707
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')
eta_squared(Tol_Chl.lme_TP3)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-------------------------------------------
Origin | 6.31e-03 | [0.00, 1.00]
Site | 0.01 | [0.00, 1.00]
Origin:Site | 2.69e-03 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
Tol_Chl_TP3.res<-data.frame(summary(Tol_Chl.lme_TP3)$coefficients[-1,])
Tol_Chl_TP3.res$Predictor<-c("Origin", "Site", "Origin x Site")
Tol_Chl_TP3.res$EtaSq<-c(eta_squared(Tol_Chl.lme_TP3)$Eta2)
Not normal, but model residuals are OK.
Effect size of Origin for each Site
##KL
Tol_Chl.lme_TP3_KL<-lmer(log(Chl.prop+1)~Origin+(1|Genotype), data=Therm_TP3[which(Therm_TP3$Site=="KL"),])
summary(Tol_Chl.lme_TP3_KL)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Chl.prop + 1) ~ Origin + (1 | Genotype)
Data: Therm_TP3[which(Therm_TP3$Site == "KL"), ]
REML criterion at convergence: -35
Scaled residuals:
Min 1Q Median 3Q Max
-2.1966 -0.5212 -0.1283 0.4729 1.9007
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.002703 0.05199
Residual 0.008473 0.09205
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.43851 0.04009 3.25876 10.938 0.00112 **
OriginTransplant 0.03181 0.03758 20.00000 0.847 0.40725
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.469
eta_squared(Tol_Chl.lme_TP3_KL)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.03 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##SS
Tol_Chl.lme_TP3_SS<-lmer(log(Chl.prop+1)~Origin+(1|Genotype), data=Therm_TP3[which(Therm_TP3$Site=="SS"),])
boundary (singular) fit: see help('isSingular')
summary(Tol_Chl.lme_TP3_SS)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Chl.prop + 1) ~ Origin + (1 | Genotype)
Data: Therm_TP3[which(Therm_TP3$Site == "SS"), ]
REML criterion at convergence: -17.1
Scaled residuals:
Min 1Q Median 3Q Max
-1.4740 -0.7462 -0.2782 0.6974 1.4761
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 3.210e-22 1.792e-11
Residual 2.144e-02 1.464e-01
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.477025 0.042265 22.000000 11.286 1.28e-10 ***
OriginTransplant 0.006709 0.059772 22.000000 0.112 0.912
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.707
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')
eta_squared(Tol_Chl.lme_TP3_SS)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 5.72e-04 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save results
Tol_Chl_TP3.site.res<-data.frame(rbind(summary(Tol_Chl.lme_TP3_KL)$coefficients[-1,],
summary(Tol_Chl.lme_TP3_SS)$coefficients[-1,]))
Tol_Chl_TP3.site.res$Predictor<-c("KL Origin", "SS Origin")
Tol_Chl_TP3.site.res$EtaSq<-c(eta_squared(Tol_Chl.lme_TP3_KL)$Eta2, eta_squared(Tol_Chl.lme_TP3_SS)$Eta2)
##Combine results
Tol_Chl_TP3.res<-rbind(Tol_Chl_TP3.res, Tol_Chl_TP3.site.res)
Tol_Chl_TP3.res$Response<-rep("Chlorophyll Retention", nrow(Tol_Chl_TP3.res))
##Summary statistics by Site and Origin
TP3_TolChl.sum<-summarySE(Therm_TP3, measurevar="Chl.prop", groupvars=c("Site", "Origin", "Site.Orig"), na.rm=TRUE)
##Plot Average Retention across Treatments
TP3_TolChl.plot<-ggplot(TP3_TolChl.sum, aes(x=Site, y=Chl.prop, colour=Site.Orig)) +
scale_colour_manual(values=Orig.colors.o)+
geom_errorbar(aes(ymin=Chl.prop-se, ymax=Chl.prop+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Origin", y="Chlorophyll Retention", colour="Origin")+
ylim(0, 1); TP3_TolChl.plot
#Check normality
hist(Therm_TP3$PAM.prop)
shapiro.test(Therm_TP3$PAM.prop)
Shapiro-Wilk normality test
data: Therm_TP3$PAM.prop
W = 0.96358, p-value = 0.1408
#Normal
##Model
#Function of Site and Origin, with Genotype as a Random effect
#Interaction between Site and Origin
Tol_PAM.lme_TP3<-lmer(PAM.prop~Origin*Site+(1|Genotype), data=Therm_TP3)
##Check residuals
Tol_PAM.lme_res_TP3 <- simulateResiduals(fittedModel = Tol_PAM.lme_TP3, plot = F)
plot(Tol_PAM.lme_res_TP3)
##Model results
summary(Tol_PAM.lme_TP3)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: PAM.prop ~ Origin * Site + (1 | Genotype)
Data: Therm_TP3
REML criterion at convergence: -172.6
Scaled residuals:
Min 1Q Median 3Q Max
-2.29399 -0.71463 0.01706 0.83761 1.77082
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 3.177e-05 0.005637
Residual 9.061e-04 0.030101
Number of obs: 48, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.954325 0.009279 14.520200 102.850 <2e-16 ***
OriginTransplant 0.009025 0.012289 42.000000 0.734 0.467
SiteSS 0.003975 0.012289 42.000000 0.323 0.748
OriginTransplant:SiteSS -0.019133 0.017379 42.000000 -1.101 0.277
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) OrgnTr SiteSS
OrgnTrnspln -0.662
SiteSS -0.662 0.500
OrgnTrn:SSS 0.468 -0.707 -0.707
eta_squared(Tol_PAM.lme_TP3)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-------------------------------------------
Origin | 9.25e-05 | [0.00, 1.00]
Site | 9.76e-03 | [0.00, 1.00]
Origin:Site | 0.03 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
Tol_PAM_TP3.res<-data.frame(summary(Tol_PAM.lme_TP3)$coefficients[-1,])
Tol_PAM_TP3.res$Predictor<-c("Origin", "Site", "Origin x Site")
Tol_PAM_TP3.res$EtaSq<-c(eta_squared(Tol_PAM.lme_TP3)$Eta2)
Effect size of Origin for each Site
##KL
Tol_PAM.lme_TP3_KL<-lmer(PAM.prop~Origin+(1|Genotype), data=Therm_TP3[which(Therm_TP3$Site=="KL"),])
summary(Tol_PAM.lme_TP3_KL)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: PAM.prop ~ Origin + (1 | Genotype)
Data: Therm_TP3[which(Therm_TP3$Site == "KL"), ]
REML criterion at convergence: -87.6
Scaled residuals:
Min 1Q Median 3Q Max
-2.4364 -0.5727 0.1072 0.8113 1.2811
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 3.996e-05 0.006322
Residual 8.452e-04 0.029072
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.954325 0.009152 5.656770 104.28 1.65e-10 ***
OriginTransplant 0.009025 0.011869 20.001145 0.76 0.456
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.648
eta_squared(Tol_PAM.lme_TP3_KL)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.03 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##SS
Tol_PAM.lme_TP3_SS<-lmer(PAM.prop~Origin+(1|Genotype), data=Therm_TP3[which(Therm_TP3$Site=="SS"),])
boundary (singular) fit: see help('isSingular')
summary(Tol_PAM.lme_TP3_SS)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: PAM.prop ~ Origin + (1 | Genotype)
Data: Therm_TP3[which(Therm_TP3$Site == "SS"), ]
REML criterion at convergence: -84.9
Scaled residuals:
Min 1Q Median 3Q Max
-1.6542 -0.7178 -0.1194 0.6304 1.6515
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.0000000 0.00000
Residual 0.0009841 0.03137
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.958300 0.009056 22.000000 105.820 <2e-16 ***
OriginTransplant -0.010108 0.012807 22.000000 -0.789 0.438
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.707
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')
eta_squared(Tol_PAM.lme_TP3_SS)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.03 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save results
Tol_PAM_TP3.site.res<-data.frame(rbind(summary(Tol_PAM.lme_TP3_KL)$coefficients[-1,],
summary(Tol_PAM.lme_TP3_SS)$coefficients[-1,]))
Tol_PAM_TP3.site.res$Predictor<-c("KL Origin", "SS Origin")
Tol_PAM_TP3.site.res$EtaSq<-c(eta_squared(Tol_PAM.lme_TP3_KL)$Eta2, eta_squared(Tol_PAM.lme_TP3_SS)$Eta2)
##Combine results
Tol_PAM_TP3.res<-rbind(Tol_PAM_TP3.res, Tol_PAM_TP3.site.res)
Tol_PAM_TP3.res$Response<-rep("FvFm Retention", nrow(Tol_PAM_TP3.res))
##Summary statistics by Site and Origin
TP3_TolPAM.sum<-summarySE(Therm_TP3, measurevar="PAM.prop", groupvars=c("Site", "Origin", "Site.Orig"), na.rm=TRUE)
##Plot Average Retention across Treatments
TP3_TolPAM.plot<-ggplot(TP3_TolPAM.sum, aes(x=Site, y=PAM.prop, colour=Site.Orig)) +
scale_colour_manual(values=Orig.colors.o)+
geom_errorbar(aes(ymin=PAM.prop-se, ymax=PAM.prop+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Origin", y="Photochemical Efficiency Retention", colour="Origin")+
ylim(0.7, 1); TP3_TolPAM.plot
##Combine Results
Tol_T3.res<-rbind(Tol_PAM_TP3.res, Tol_Sym_TP3.res, Tol_Chl_TP3.res)
##Add Timepoint
Tol_T3.res$TimeP<-rep("TP3", nrow(Tol_T3.res))
##Subset Timepoint 4
Therm_TP4<-subset(Therm_H, TimeP=="TP4")
#Check normality
hist(Therm_TP4$Sym.prop)
shapiro.test(Therm_TP4$Sym.prop)
Shapiro-Wilk normality test
data: Therm_TP4$Sym.prop
W = 0.93453, p-value = 0.01005
#Not normal
hist(log(Therm_TP4$Sym.prop+1))
shapiro.test(log(Therm_TP4$Sym.prop+1))
Shapiro-Wilk normality test
data: log(Therm_TP4$Sym.prop + 1)
W = 0.95458, p-value = 0.061
#Normal
##Model with log +1 transformation
#Function of Site and Origin, with Genotype as a Random effect
#Interaction between Site and Origin
Tol_Sym.lme_TP4<-lmer(log(Sym.prop+1)~Origin*Site+(1|Genotype), data=Therm_TP4)
##Check residuals
Tol_Sym.lme_res_TP4 <- simulateResiduals(fittedModel = Tol_Sym.lme_TP4, plot = F)
plot(Tol_Sym.lme_res_TP4)
##Model results
summary(Tol_Sym.lme_TP4)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Sym.prop + 1) ~ Origin * Site + (1 | Genotype)
Data: Therm_TP4
REML criterion at convergence: -44.7
Scaled residuals:
Min 1Q Median 3Q Max
-1.5412 -0.7855 -0.1534 0.5518 2.2412
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.001295 0.03598
Residual 0.016289 0.12763
Number of obs: 48, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.24607 0.04230 9.94436 5.818 0.000173 ***
OriginTransplant 0.04098 0.05210 42.00002 0.787 0.435975
SiteSS -0.03092 0.05210 42.00002 -0.593 0.556032
OriginTransplant:SiteSS -0.03256 0.07369 42.00002 -0.442 0.660816
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) OrgnTr SiteSS
OrgnTrnspln -0.616
SiteSS -0.616 0.500
OrgnTrn:SSS 0.436 -0.707 -0.707
eta_squared(Tol_Sym.lme_TP4)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-------------------------------------------
Origin | 0.01 | [0.00, 1.00]
Site | 0.04 | [0.00, 1.00]
Origin:Site | 4.63e-03 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
Tol_Sym_TP4.res<-data.frame(summary(Tol_Sym.lme_TP4)$coefficients[-1,])
Tol_Sym_TP4.res$Predictor<-c("Origin", "Site", "Origin x Site")
Tol_Sym_TP4.res$EtaSq<-c(eta_squared(Tol_Sym.lme_TP4)$Eta2)
Effect size of Origin for each Site
##KL
Tol_Sym.lme_TP4_KL<-lmer(log(Sym.prop+1)~Origin+(1|Genotype), data=Therm_TP4[which(Therm_TP4$Site=="KL"),])
summary(Tol_Sym.lme_TP4_KL)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Sym.prop + 1) ~ Origin + (1 | Genotype)
Data: Therm_TP4[which(Therm_TP4$Site == "KL"), ]
REML criterion at convergence: -26.2
Scaled residuals:
Min 1Q Median 3Q Max
-1.3745 -0.8456 0.0271 0.5432 2.0020
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.002906 0.05391
Residual 0.012909 0.11362
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.24607 0.04522 3.63655 5.442 0.00725 **
OriginTransplant 0.04098 0.04638 20.00000 0.884 0.38746
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.513
eta_squared(Tol_Sym.lme_TP4_KL)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.04 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##SS
Tol_Sym.lme_TP4_SS<-lmer(log(Sym.prop+1)~Origin+(1|Genotype), data=Therm_TP4[which(Therm_TP4$Site=="SS"),])
boundary (singular) fit: see help('isSingular')
summary(Tol_Sym.lme_TP4_SS)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Sym.prop + 1) ~ Origin + (1 | Genotype)
Data: Therm_TP4[which(Therm_TP4$Site == "SS"), ]
REML criterion at convergence: -19.3
Scaled residuals:
Min 1Q Median 3Q Max
-1.19030 -0.85642 -0.01412 0.53810 2.27217
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.00000 0.0000
Residual 0.01944 0.1394
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.215142 0.040247 22.000000 5.346 2.29e-05 ***
OriginTransplant 0.008418 0.056918 22.000000 0.148 0.884
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.707
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')
eta_squared(Tol_Sym.lme_TP4_SS)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 9.93e-04 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save results
Tol_Sym_TP4.site.res<-data.frame(rbind(summary(Tol_Sym.lme_TP4_KL)$coefficients[-1,],
summary(Tol_Sym.lme_TP4_SS)$coefficients[-1,]))
Tol_Sym_TP4.site.res$Predictor<-c("KL Origin", "SS Origin")
Tol_Sym_TP4.site.res$EtaSq<-c(eta_squared(Tol_Sym.lme_TP4_KL)$Eta2, eta_squared(Tol_Sym.lme_TP4_SS)$Eta2)
##Combine results
Tol_Sym_TP4.res<-rbind(Tol_Sym_TP4.res, Tol_Sym_TP4.site.res)
Tol_Sym_TP4.res$Response<-rep("Symbiont Retention", nrow(Tol_Sym_TP4.res))
##Summary statistics by Site and Origin
TP4_TolSym.sum<-summarySE(Therm_TP4, measurevar="Sym.prop", groupvars=c("Site", "Origin", "Site.Orig"), na.rm=TRUE)
##Plot Average Retention across Treatments
TP4_TolSym.plot<-ggplot(TP4_TolSym.sum, aes(x=Site, y=Sym.prop, colour=Site.Orig)) +
scale_colour_manual(values=Orig.colors.o)+
geom_errorbar(aes(ymin=Sym.prop-se, ymax=Sym.prop+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Origin", y="Symbiont Retention", colour="Origin")+
ylim(0, 1); TP4_TolSym.plot
#Check normality
hist(Therm_TP4$Chl.prop)
shapiro.test(Therm_TP4$Chl.prop)
Shapiro-Wilk normality test
data: Therm_TP4$Chl.prop
W = 0.94204, p-value = 0.02131
#Not normal
hist(log(Therm_TP4$Chl.prop+1))
shapiro.test(log(Therm_TP4$Chl.prop+1))
Shapiro-Wilk normality test
data: log(Therm_TP4$Chl.prop + 1)
W = 0.9512, p-value = 0.04824
##Nearly normal
##Model with log +1 transformation
#Function of Site and Origin, with Genotype as a Random effect
#Interaction between Site and Origin
Tol_Chl.lme_TP4<-lmer(log(Chl.prop+1)~Origin*Site+(1|Genotype), data=Therm_TP4)
boundary (singular) fit: see help('isSingular')
##Check residuals
Tol_Chl.lme_res_TP4 <- simulateResiduals(fittedModel = Tol_Chl.lme_TP4, plot = F)
plot(Tol_Chl.lme_res_TP4)
##Model results
summary(Tol_Chl.lme_TP4)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Chl.prop + 1) ~ Origin * Site + (1 | Genotype)
Data: Therm_TP4
REML criterion at convergence: -139.8
Scaled residuals:
Min 1Q Median 3Q Max
-1.6397 -0.7647 -0.1492 0.5006 2.4251
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.000000 0.00000
Residual 0.001803 0.04246
Number of obs: 47, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.09448 0.01226 43.00000 7.707 1.24e-09 ***
OriginTransplant -0.01266 0.01734 43.00000 -0.730 0.469
SiteSS -0.02659 0.01773 43.00000 -1.500 0.141
OriginTransplant:SiteSS 0.02939 0.02479 43.00000 1.185 0.242
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) OrgnTr SiteSS
OrgnTrnspln -0.707
SiteSS -0.692 0.489
OrgnTrn:SSS 0.494 -0.699 -0.715
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')
eta_squared(Tol_Chl.lme_TP4)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-------------------------------------------
Origin | 6.26e-04 | [0.00, 1.00]
Site | 0.02 | [0.00, 1.00]
Origin:Site | 0.03 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
Tol_Chl_TP4.res<-data.frame(summary(Tol_Chl.lme_TP4)$coefficients[-1,])
Tol_Chl_TP4.res$Predictor<-c("Origin", "Site", "Origin x Site")
Tol_Chl_TP4.res$EtaSq<-c(eta_squared(Tol_Chl.lme_TP4)$Eta2)
Effect size of Origin for each Site
##KL
Tol_Chl.lme_TP4_KL<-lmer(log(Chl.prop+1)~Origin+(1|Genotype), data=Therm_TP4[which(Therm_TP4$Site=="KL"),])
boundary (singular) fit: see help('isSingular')
summary(Tol_Chl.lme_TP4_KL)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Chl.prop + 1) ~ Origin + (1 | Genotype)
Data: Therm_TP4[which(Therm_TP4$Site == "KL"), ]
REML criterion at convergence: -76.7
Scaled residuals:
Min 1Q Median 3Q Max
-1.5133 -0.5558 -0.1258 0.5428 2.7208
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.000000 0.00000
Residual 0.001433 0.03785
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.09448 0.01093 22.00000 8.647 1.59e-08 ***
OriginTransplant -0.01266 0.01545 22.00000 -0.819 0.421
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.707
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')
eta_squared(Tol_Chl.lme_TP4_KL)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.03 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##SS
Tol_Chl.lme_TP4_SS<-lmer(log(Chl.prop+1)~Origin+(1|Genotype), data=Therm_TP4[which(Therm_TP4$Site=="SS"),])
boundary (singular) fit: see help('isSingular')
summary(Tol_Chl.lme_TP4_SS)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: log(Chl.prop + 1) ~ Origin + (1 | Genotype)
Data: Therm_TP4[which(Therm_TP4$Site == "SS"), ]
REML criterion at convergence: -64.1
Scaled residuals:
Min 1Q Median 3Q Max
-1.4873 -0.7211 -0.1933 0.5918 2.1830
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.000000 0.00000
Residual 0.002192 0.04681
Number of obs: 23, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.06788 0.01412 21.00000 4.809 9.41e-05 ***
OriginTransplant 0.01673 0.01954 21.00000 0.856 0.402
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.722
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')
eta_squared(Tol_Chl.lme_TP4_SS)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.03 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save results
Tol_Chl_TP4.site.res<-data.frame(rbind(summary(Tol_Chl.lme_TP4_KL)$coefficients[-1,],
summary(Tol_Chl.lme_TP4_SS)$coefficients[-1,]))
Tol_Chl_TP4.site.res$Predictor<-c("KL Origin", "SS Origin")
Tol_Chl_TP4.site.res$EtaSq<-c(eta_squared(Tol_Chl.lme_TP4_KL)$Eta2, eta_squared(Tol_Chl.lme_TP4_SS)$Eta2)
##Combine results
Tol_Chl_TP4.res<-rbind(Tol_Chl_TP4.res, Tol_Chl_TP4.site.res)
Tol_Chl_TP4.res$Response<-rep("Chlorophyll Retention", nrow(Tol_Chl_TP4.res))
##Summary statistics by Site and Origin
TP4_TolChl.sum<-summarySE(Therm_TP4, measurevar="Chl.prop", groupvars=c("Site", "Origin", "Site.Orig"), na.rm=TRUE)
##Plot Average Retention across Treatments
TP4_TolChl.plot<-ggplot(TP4_TolChl.sum, aes(x=Site, y=Chl.prop, colour=Site.Orig)) +
scale_colour_manual(values=Orig.colors.o)+
geom_errorbar(aes(ymin=Chl.prop-se, ymax=Chl.prop+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Origin", y="Chlorophyll Retention", colour="Origin")+
ylim(0, 0.35); TP4_TolChl.plot
#Check normality
hist(Therm_TP4$PAM.prop)
shapiro.test(Therm_TP4$PAM.prop)
Shapiro-Wilk normality test
data: Therm_TP4$PAM.prop
W = 0.86379, p-value = 5.159e-05
#Not normal
hist(log(Therm_TP4$PAM.prop+1))
shapiro.test(log(Therm_TP4$PAM.prop+1))
Shapiro-Wilk normality test
data: log(Therm_TP4$PAM.prop + 1)
W = 0.84113, p-value = 1.275e-05
##Still not normal
##Model
#Function of Site and Origin, with Genotype as a Random effect
#Interaction between Site and Origin
Tol_PAM.lme_TP4<-lmer(PAM.prop~Origin*Site+(1|Genotype), data=Therm_TP4)
##Check residuals
Tol_PAM.lme_res_TP4 <- simulateResiduals(fittedModel = Tol_PAM.lme_TP4, plot = F)
plot(Tol_PAM.lme_res_TP4)
##Model results
summary(Tol_PAM.lme_TP4)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: PAM.prop ~ Origin * Site + (1 | Genotype)
Data: Therm_TP4
REML criterion at convergence: -60.9
Scaled residuals:
Min 1Q Median 3Q Max
-2.1005 -0.4066 -0.0093 0.5802 2.4408
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.01070 0.1035
Residual 0.01028 0.1014
Number of obs: 48, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.74804 0.06651 2.73320 11.247 0.00228 **
OriginTransplant 0.03239 0.04138 42.00000 0.783 0.43818
SiteSS 0.02783 0.04138 42.00000 0.672 0.50502
OriginTransplant:SiteSS -0.05803 0.05852 42.00000 -0.991 0.32713
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) OrgnTr SiteSS
OrgnTrnspln -0.311
SiteSS -0.311 0.500
OrgnTrn:SSS 0.220 -0.707 -0.707
eta_squared(Tol_PAM.lme_TP4)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-------------------------------------------
Origin | 3.17e-04 | [0.00, 1.00]
Site | 3.92e-05 | [0.00, 1.00]
Origin:Site | 0.02 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save model results
Tol_PAM_TP4.res<-data.frame(summary(Tol_PAM.lme_TP4)$coefficients[-1,])
Tol_PAM_TP4.res$Predictor<-c("Origin", "Site", "Origin x Site")
Tol_PAM_TP4.res$EtaSq<-c(eta_squared(Tol_PAM.lme_TP4)$Eta2)
Effect size of Origin for each Site
##KL
Tol_PAM.lme_TP4_KL<-lmer(PAM.prop~Origin+(1|Genotype), data=Therm_TP4[which(Therm_TP4$Site=="KL"),])
summary(Tol_PAM.lme_TP4_KL)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: PAM.prop ~ Origin + (1 | Genotype)
Data: Therm_TP4[which(Therm_TP4$Site == "KL"), ]
REML criterion at convergence: -30.9
Scaled residuals:
Min 1Q Median 3Q Max
-1.84954 -0.38862 -0.07926 0.54042 2.60201
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.011716 0.10824
Residual 0.009227 0.09606
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.74804 0.06837 2.37263 10.941 0.00438 **
OriginTransplant 0.03239 0.03922 20.00001 0.826 0.41856
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.287
eta_squared(Tol_PAM.lme_TP4_KL)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.03 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##SS
Tol_PAM.lme_TP4_SS<-lmer(PAM.prop~Origin+(1|Genotype), data=Therm_TP4[which(Therm_TP4$Site=="SS"),])
summary(Tol_PAM.lme_TP4_SS)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: PAM.prop ~ Origin + (1 | Genotype)
Data: Therm_TP4[which(Therm_TP4$Site == "SS"), ]
REML criterion at convergence: -25.7
Scaled residuals:
Min 1Q Median 3Q Max
-2.08840 -0.38585 0.06028 0.59276 2.02689
Random effects:
Groups Name Variance Std.Dev.
Genotype (Intercept) 0.008458 0.09197
Residual 0.012219 0.11054
Number of obs: 24, groups: Genotype, 3
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.77587 0.06195 2.65242 12.524 0.00196 **
OriginTransplant -0.02563 0.04513 20.00000 -0.568 0.57634
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
OrgnTrnspln -0.364
eta_squared(Tol_PAM.lme_TP4_SS)
# Effect Size for ANOVA (Type III)
Parameter | Eta2 (partial) | 95% CI
-----------------------------------------
Origin | 0.02 | [0.00, 1.00]
- One-sided CIs: upper bound fixed at [1.00].
##Save results
Tol_PAM_TP4.site.res<-data.frame(rbind(summary(Tol_PAM.lme_TP4_KL)$coefficients[-1,],
summary(Tol_PAM.lme_TP4_SS)$coefficients[-1,]))
Tol_PAM_TP4.site.res$Predictor<-c("KL Origin", "SS Origin")
Tol_PAM_TP4.site.res$EtaSq<-c(eta_squared(Tol_PAM.lme_TP4_KL)$Eta2, eta_squared(Tol_PAM.lme_TP4_SS)$Eta2)
##Combine results
Tol_PAM_TP4.res<-rbind(Tol_PAM_TP4.res, Tol_PAM_TP4.site.res)
Tol_PAM_TP4.res$Response<-rep("FvFm Retention", nrow(Tol_PAM_TP4.res))
##Summary statistics by Site and Origin
TP4_TolPAM.sum<-summarySE(Therm_TP4, measurevar="PAM.prop", groupvars=c("Site", "Origin", "Site.Orig"), na.rm=TRUE)
##Plot Average Retention across Treatments
TP4_TolPAM.plot<-ggplot(TP4_TolPAM.sum, aes(x=Site, y=PAM.prop, colour=Site.Orig)) +
scale_colour_manual(values=Orig.colors.o)+
geom_errorbar(aes(ymin=PAM.prop-se, ymax=PAM.prop+se), width=cap.sz, linewidth=bar.sz, position=position_dodge(width=0.5)) +
geom_point(size=point.sz, position=position_dodge(width=0.5))+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Site and Origin", y="Photochemical Efficiency Retention", colour="Origin")+
ylim(0.6, 1); TP4_TolPAM.plot
##Combine Results
Tol_T4.res<-rbind(Tol_PAM_TP4.res, Tol_Sym_TP4.res, Tol_Chl_TP4.res)
##Add Timepoint
Tol_T4.res$TimeP<-rep("TP4", nrow(Tol_T4.res))
##Dataframe of effect size results
Tol.ES<-rbind(Tol_T1.res, Tol_T2.res, Tol_T2.res, Tol_T4.res)
Tol.ES<-Tol.ES %>% dplyr::rename(p = Pr...t..)
Tol.ES<-Tol.ES[,c("TimeP", "Predictor", "Response", "EtaSq", "p")]
##Site specific results
Tol.ES<-Tol.ES[which(Tol.ES$Predictor=="KL Origin" | Tol.ES$Predictor=="SS Origin"),]
Tol.ES<-Tol.ES %>% separate_wider_delim(cols="Predictor", delim=" ", names=c("Site", "Predictor"), cols_remove = TRUE)
##Add Metric names
Tol.ES$Metric<-ifelse(Tol.ES$Response== "FvFm Retention", "FvFm", ifelse(
Tol.ES$Response== "Symbiont Retention", "Symbionts", ifelse(
Tol.ES$Response=="Chlorophyll Retention", "Chlorophyll", NA)))
##Add Significance levels
Tol.ES$Sig<-ifelse(Tol.ES$p<0.001, "***", ifelse(Tol.ES$p<0.01, "**",
ifelse(Tol.ES$p<0.05, "*", ifelse(Tol.ES$p<0.1, "-", NA))))
Tol_Sym.ES.plot<-ggplot(Tol.ES[which(Tol.ES$Metric=="Symbionts"),], aes(x=TimeP, y=EtaSq, fill=Site))+
geom_bar(stat="identity", position=position_dodge())+
scale_fill_manual(values=Site.colors.o)+
theme_classic()+
ggtitle("Symbiont Retention")+
theme(plot.title = element_text(colour="black", size=panel.lab.sz, hjust=0.5), axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"), legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Time Point", y=expression(paste("Effect Size (p", eta^2, ")")), colour="Site")+
ylim(0, 0.4)+
geom_text(aes(label=Sig), vjust=-0.02, color="black", position=position_dodge(0.9), size=levels.sz, fontface="bold"); Tol_Sym.ES.plot
Tol_Chl.ES.plot<-ggplot(Tol.ES[which(Tol.ES$Metric=="Chlorophyll"),], aes(x=TimeP, y=EtaSq, fill=Site))+
geom_bar(stat="identity", position=position_dodge())+
scale_fill_manual(values=Site.colors.o)+
theme_classic()+
ggtitle("Chlorophyll Retention")+
theme(plot.title = element_text(colour="black", size=panel.lab.sz, hjust=0.5), axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"), legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Time Point", y=expression(paste("Effect Size (p", eta^2, ")")), colour="Site")+
ylim(0, 0.4)+
geom_text(aes(label=Sig), vjust=-0.02, color="black", position=position_dodge(0.9), size=levels.sz, fontface="bold"); Tol_Chl.ES.plot
Tol_PAM.ES.plot<-ggplot(Tol.ES[which(Tol.ES$Metric=="FvFm"),], aes(x=TimeP, y=EtaSq, fill=Site))+
geom_bar(stat="identity", position=position_dodge())+
scale_fill_manual(values=Site.colors.o)+
theme_classic()+
ggtitle("Photochemical Efficiency Retention")+
theme(plot.title = element_text(colour="black", size=panel.lab.sz, hjust=0.5), axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"), legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x="Time Point", y=expression(paste("Effect Size (p", eta^2, ")")), colour="Site")+
ylim(0, 0.4)+
geom_text(aes(label=Sig), vjust=-0.02, color="black", position=position_dodge(0.9), size=levels.sz, fontface="bold"); Tol_PAM.ES.plot
Calculating the average of Growth (TP1 to TP2) and Thermal Tolerance (Retention of Symbionts, Chlorophyll, and Photochemical Efficiency at TP1) for each Set (Site, Genotype, Origin).
##Average
names(Growth_T1.2)
[1] "Site" "ID" "Orig" "Origin" "Genotype" "TL.TP1_cm" "TL.TP2_cm"
[8] "TL.TP3_cm" "TL.TP4_cm" "Set" "Site.Orig" "TP1" "TP2" "TP1.2_days"
[15] "Ext_cm" "TLE_cm.day"
T1.2_Growth_Set<-aggregate(Growth_T1.2$TLE_cm.day, list(Growth_T1.2$Site, Growth_T1.2$Genotype, Growth_T1.2$Orig), mean)
names(T1.2_Growth_Set)<-c("Site", "Genotype", "Orig", "TLE_cm.day")
##Add Standard Error
T1.2_Growth_Set$TLE_se<-aggregate(Growth_T1.2$TLE_cm.day, list(Growth_T1.2$Site, Growth_T1.2$Genotype, Growth_T1.2$Orig), std.error)[,4]
##Average
names(Therm_TP1)
[1] "TimeP" "Site" "Genotype" "Origin" "ID" "RandN"
[7] "Treat" "Treatment" "Orig" "Set" "Site.Orig" "SA_cm2"
[13] "Chl_ug.cm2" "Sym10.6_cm2" "Fv_Fm" "Chl_ug.cm2_C" "Sym10.6_cm2_C" "Fv_Fm_C"
[19] "Chl.prop" "Sym.prop" "PAM.prop"
T1_Therm_Set<-aggregate(Therm_TP1[,c(19:21)], list(Therm_TP1$Site, Therm_TP1$Genotype, Therm_TP1$Orig), mean)
names(T1_Therm_Set)<-c("Site", "Genotype", "Orig", "Chl.prop", "Sym.prop", "PAM.prop")
##Add Standard Error
T1_Therm_Set$Chl_se<-aggregate(Therm_TP1$Chl.prop, list(Therm_TP1$Site, Therm_TP1$Genotype, Therm_TP1$Orig), std.error)[,4]
T1_Therm_Set$Sym_se<-aggregate(Therm_TP1$Sym.prop, list(Therm_TP1$Site, Therm_TP1$Genotype, Therm_TP1$Orig), std.error)[,4]
T1_Therm_Set$PAM_se<-aggregate(Therm_TP1$PAM.prop, list(Therm_TP1$Site, Therm_TP1$Genotype, Therm_TP1$Orig), std.error)[,4]
T1.2_Perf_Set<-merge(T1.2_Growth_Set, T1_Therm_Set, all=TRUE)
T1.2_Perf_Set$Set<-paste(T1.2_Perf_Set$Site, T1.2_Perf_Set$Genotype, T1.2_Perf_Set$Orig, sep=".")
##Set factor variables
T1.2_Perf_Set$Site<-factor(T1.2_Perf_Set$Site, levels=c("KL", "SS"))
T1.2_Perf_Set$Genotype<-factor(T1.2_Perf_Set$Genotype, levels=c("AC8", "AC10", "AC12"))
T1.2_Perf_Set$Orig<-factor(T1.2_Perf_Set$Orig, levels=c("N", "T"))
##Add a Sample Set Variable
T1.2_Perf_Set$Set<-paste(T1.2_Perf_Set$Site, T1.2_Perf_Set$Genotype, T1.2_Perf_Set$Orig, sep=".")
T1.2_Perf_Set$Set<-factor(T1.2_Perf_Set$Set)
##Add Site.Orig variable
T1.2_Perf_Set$Site.Orig<-paste(T1.2_Perf_Set$Site, T1.2_Perf_Set$Orig, sep=".")
T1.2_Perf_Set$Site.Orig<-factor(T1.2_Perf_Set$Site.Orig, levels=c("KL.N", "KL.T", "SS.N", "SS.T"))
cor.test(T1.2_Perf_Set$TLE_cm.day, T1.2_Perf_Set$Sym.prop, method="spearman")
Spearman's rank correlation rho
data: T1.2_Perf_Set$TLE_cm.day and T1.2_Perf_Set$Sym.prop
S = 462, p-value = 0.03733
alternative hypothesis: true rho is not equal to 0
sample estimates:
rho
-0.6153846
T1.2_TLE_Sym.plot<-ggplot(T1.2_Perf_Set, aes(x=TLE_cm.day, y=Sym.prop)) +
geom_smooth(method="lm", color="#AA185AFF", fill="#F8D7BFFF", se=TRUE)+
geom_point(aes(colour=Site.Orig, shape=Genotype),size=point.sz)+
scale_colour_manual(values=Orig.colors.o)+
scale_shape_manual(values=Geno.shapes.o)+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x=expression(paste('Growth Rate (cm day'^-1*")")),
y="Symbiont Retention",
colour="Site Origin")+
annotate("text", x = 0.07, y = -.15, label=expression(bolditalic(paste(r[S], " = -0.615, p = 0.037"))), size=sig.sz, hjust = 0); T1.2_TLE_Sym.plot
cor.test(T1.2_Perf_Set$TLE_cm.day, T1.2_Perf_Set$Chl.prop, method="spearman")
Spearman's rank correlation rho
data: T1.2_Perf_Set$TLE_cm.day and T1.2_Perf_Set$Chl.prop
S = 454, p-value = 0.04884
alternative hypothesis: true rho is not equal to 0
sample estimates:
rho
-0.5874126
T1.2_TLE_Chl.plot<-ggplot(T1.2_Perf_Set, aes(x=TLE_cm.day, y=Chl.prop)) +
geom_smooth(method="lm", color="#AA185AFF", fill="#F8D7BFFF", se=TRUE)+
geom_point(aes(colour=Site.Orig, shape=Genotype),size=point.sz)+
scale_colour_manual(values=Orig.colors.o)+
scale_shape_manual(values=Geno.shapes.o)+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x=expression(paste('Growth Rate (cm day'^-1*")")),
y="Chlorophyll Retention",
colour="Site Origin")+
annotate("text", x = 0.07, y = -.05, label=expression(bolditalic(paste(r[S], " = -0.587, p = 0.049"))), size=sig.sz, hjust = 0); T1.2_TLE_Chl.plot
cor.test(T1.2_Perf_Set$TLE_cm.day, T1.2_Perf_Set$PAM.prop, method="spearman")
Spearman's rank correlation rho
data: T1.2_Perf_Set$TLE_cm.day and T1.2_Perf_Set$PAM.prop
S = 480, p-value = 0.01883
alternative hypothesis: true rho is not equal to 0
sample estimates:
rho
-0.6783217
T1.2_TLE_PAM.plot<-ggplot(T1.2_Perf_Set, aes(x=TLE_cm.day, y=PAM.prop)) +
geom_smooth(method="lm", color="#AA185AFF", fill="#F8D7BFFF", se=TRUE)+
geom_point(aes(colour=Site.Orig, shape=Genotype),size=point.sz)+
scale_colour_manual(values=Orig.colors.o)+
scale_shape_manual(values=Geno.shapes.o)+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x=expression(paste('Growth Rate (cm day'^-1*")")),
y="Photo. Effic. Retention",
colour="Site Origin")+
annotate("text", x = 0.07, y = 0.7, label=expression(bolditalic(paste(r[S], " = -0.678, p = 0.019"))), size=sig.sz, hjust = 0); T1.2_TLE_PAM.plot
Calculating the average of Growth (TP3 to TP4) and Thermal Tolerance (Retention of Symbionts, Chlorophyll, and Photochemical Efficiency at TP3) for each Set (Site, Genotype, Origin).
##Average
names(Growth_T3.4)
[1] "Site" "ID" "Orig" "Origin" "Genotype" "TL.TP1_cm" "TL.TP2_cm"
[8] "TL.TP3_cm" "TL.TP4_cm" "Set" "Site.Orig" "TP3" "TP4" "TP3.4_days"
[15] "Ext_cm" "TLE_cm.day"
T3.4_Growth_Set<-aggregate(Growth_T3.4$TLE_cm.day, list(Growth_T3.4$Site, Growth_T3.4$Genotype, Growth_T3.4$Orig), mean)
names(T3.4_Growth_Set)<-c("Site", "Genotype", "Orig", "TLE_cm.day")
##Add Standard Error
T3.4_Growth_Set$TLE_se<-aggregate(Growth_T3.4$TLE_cm.day, list(Growth_T3.4$Site, Growth_T3.4$Genotype, Growth_T3.4$Orig), std.error)[,4]
##Average
names(Therm_TP4)
[1] "TimeP" "Site" "Genotype" "Origin" "ID" "RandN"
[7] "Treat" "Treatment" "Orig" "Set" "Site.Orig" "SA_cm2"
[13] "Chl_ug.cm2" "Sym10.6_cm2" "Fv_Fm" "Chl_ug.cm2_C" "Sym10.6_cm2_C" "Fv_Fm_C"
[19] "Chl.prop" "Sym.prop" "PAM.prop"
T4_Therm_Set<-aggregate(Therm_TP4[,c(19:21)], list(Therm_TP4$Site, Therm_TP4$Genotype, Therm_TP4$Orig), mean)
names(T4_Therm_Set)<-c("Site", "Genotype", "Orig", "Chl.prop", "Sym.prop", "PAM.prop")
##Add Standard Error
T4_Therm_Set$Chl_se<-aggregate(Therm_TP4$Chl.prop, list(Therm_TP4$Site, Therm_TP4$Genotype, Therm_TP4$Orig), std.error)[,4]
T4_Therm_Set$Sym_se<-aggregate(Therm_TP4$Sym.prop, list(Therm_TP4$Site, Therm_TP4$Genotype, Therm_TP4$Orig), std.error)[,4]
T4_Therm_Set$PAM_se<-aggregate(Therm_TP4$PAM.prop, list(Therm_TP4$Site, Therm_TP4$Genotype, Therm_TP4$Orig), std.error)[,4]
T3.4_Perf_Set<-merge(T3.4_Growth_Set, T4_Therm_Set, all=TRUE)
T3.4_Perf_Set$Set<-paste(T3.4_Perf_Set$Site, T3.4_Perf_Set$Genotype, T3.4_Perf_Set$Orig, sep=".")
##Set factor variables
T3.4_Perf_Set$Site<-factor(T3.4_Perf_Set$Site, levels=c("KL", "SS"))
T3.4_Perf_Set$Genotype<-factor(T3.4_Perf_Set$Genotype, levels=c("AC8", "AC10", "AC12"))
T3.4_Perf_Set$Orig<-factor(T3.4_Perf_Set$Orig, levels=c("N", "T"))
##Add a Sample Set Variable
T3.4_Perf_Set$Set<-paste(T3.4_Perf_Set$Site, T3.4_Perf_Set$Genotype, T3.4_Perf_Set$Orig, sep=".")
T3.4_Perf_Set$Set<-factor(T3.4_Perf_Set$Set)
##Add Site.Orig variable
T3.4_Perf_Set$Site.Orig<-paste(T3.4_Perf_Set$Site, T3.4_Perf_Set$Orig, sep=".")
T3.4_Perf_Set$Site.Orig<-factor(T3.4_Perf_Set$Site.Orig, levels=c("KL.N", "KL.T", "SS.N", "SS.T"))
cor.test(T3.4_Perf_Set$TLE_cm.day, T3.4_Perf_Set$Sym.prop, method="spearman")
Spearman's rank correlation rho
data: T3.4_Perf_Set$TLE_cm.day and T3.4_Perf_Set$Sym.prop
S = 444, p-value = 0.06663
alternative hypothesis: true rho is not equal to 0
sample estimates:
rho
-0.5524476
T3.4_TLE_Sym.plot<-ggplot(T3.4_Perf_Set, aes(x=TLE_cm.day, y=Sym.prop)) +
geom_smooth(method="lm", color="#AA185AFF", fill="#F8D7BFFF", se=TRUE)+
geom_point(aes(colour=Site.Orig, shape=Genotype),size=point.sz)+
scale_colour_manual(values=Orig.colors.o)+
scale_shape_manual(values=Geno.shapes.o)+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x=expression(paste('Growth Rate (cm day'^-1*")")),
y="Symbiont Retention",
colour="Site Origin")+
annotate("text", x = 0.07, y = 0.1, label=expression(italic(paste(r[S], " = -0.552, p = 0.067"))), size=sig.sz, hjust = 0); T3.4_TLE_Sym.plot
cor.test(T3.4_Perf_Set$TLE_cm.day, T3.4_Perf_Set$Chl.prop, method="spearman")
Spearman's rank correlation rho
data: T3.4_Perf_Set$TLE_cm.day and T3.4_Perf_Set$Chl.prop
S = 148, p-value = 0.327
alternative hypothesis: true rho is not equal to 0
sample estimates:
rho
0.3272727
T3.4_TLE_Chl.plot<-ggplot(T3.4_Perf_Set, aes(x=TLE_cm.day, y=Chl.prop)) +
geom_smooth(method="lm", color="#AA185AFF", fill="#F8D7BFFF", se=TRUE)+
geom_point(aes(colour=Site.Orig, shape=Genotype),size=point.sz)+
scale_colour_manual(values=Orig.colors.o)+
scale_shape_manual(values=Geno.shapes.o)+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x=expression(paste('Growth Rate (cm day'^-1*")")),
y="Chlorophyll Retention",
colour="Site Origin")+
annotate("text", x = 0.07, y = 0.065, label=expression(italic(paste(r[S], " = 0.327, p = 0.327"))), size=sig.sz, hjust = 0); T3.4_TLE_Chl.plot
cor.test(T3.4_Perf_Set$TLE_cm.day, T3.4_Perf_Set$PAM.prop, method="spearman")
Spearman's rank correlation rho
data: T3.4_Perf_Set$TLE_cm.day and T3.4_Perf_Set$PAM.prop
S = 500, p-value = 0.007353
alternative hypothesis: true rho is not equal to 0
sample estimates:
rho
-0.7482517
T3.4_TLE_PAM.plot<-ggplot(T3.4_Perf_Set, aes(x=TLE_cm.day, y=PAM.prop)) +
geom_smooth(method="lm", color="#AA185AFF", fill="#F8D7BFFF", se=TRUE)+
geom_point(aes(colour=Site.Orig, shape=Genotype),size=point.sz)+
scale_colour_manual(values=Orig.colors.o)+
scale_shape_manual(values=Geno.shapes.o)+
theme_classic()+
theme(axis.title.x = element_text(size = axis.title.sz), axis.title.y = element_text(size = axis.title.sz),
axis.text.x=element_text(size=axis.txt.sz, colour="black"), axis.text.y=element_text(size=axis.txt.sz, colour="black"),
legend.text=element_text(size=leg.txt.sz), legend.title=element_text(size=leg.title.sz),
legend.box.background = element_rect(color = "black"))+
labs(x=expression(paste('Growth Rate (cm day'^-1*")")),
y="Photo. Effic. Retention",
colour="Site Origin")+
annotate("text", x = 0.07, y = 0.5, label=expression(bolditalic(paste(r[S], " = -0.748, p = 0.007"))), size=sig.sz, hjust = 0); T3.4_TLE_PAM.plot
##TP 1 to TP 2
TP1.2_Growth.plot<-TP1.2_Growth.plot+ ggtitle("TP1 to TP2")+
theme(plot.title = element_text(colour="black",
size=panel.lab.sz, hjust=0.5),
legend.position=c(0.25, 0.85), legend.direction="horizontal")+
labs(colour=NULL)+
guides(color=guide_legend(nrow=2, byrow=FALSE))
##TP 3 to TP 4
TP3.4_Growth.plot<-TP3.4_Growth.plot+ ggtitle("TP3 to TP4")+
theme(plot.title = element_text(colour="black",
size=panel.lab.sz, hjust=0.5),
legend.position="none")
##Create Panel
Growth_fig<-plot_grid(TP1.2_Growth.plot, TP3.4_Growth.plot,
rel_widths=1, rel_heights = 1,
nrow=1, ncol=2, byrow=T, labels = c("A", "B"),
label_size=panel.lab.sz)
##Save Figure
ggsave(filename="Figures/03_Performance/Fig4_Growth.png", plot=Growth_fig, dpi=300, width=10, height=4.5, units="in")
##FvFm Retention
TP1_TolPAM.plot<-TP1_TolPAM.plot + ggtitle("TP1")+
labs(x="", y="Photo. Effic. Retention", colour=NULL)+
theme(plot.title = element_text(colour="black", size=panel.lab.sz, face="bold", hjust=0.5),
legend.position=c(0.5, 0.15), legend.direction="horizontal")+
guides(color=guide_legend(nrow=2, byrow=FALSE))
TP2_TolPAM.plot<-TP2_TolPAM.plot+ ggtitle("TP2")+ labs(x="", y="") +
theme(plot.title = element_text(colour="black", size=panel.lab.sz, face="bold", hjust=0.5),
legend.position="none")
TP3_TolPAM.plot<-TP3_TolPAM.plot+ ggtitle("TP3")+ labs(x="", y="") +
theme(plot.title = element_text(colour="black", size=panel.lab.sz, face="bold", hjust=0.5),
legend.position="none")
TP4_TolPAM.plot<-TP4_TolPAM.plot+ ggtitle("TP4")+ labs(x="", y="") +
theme(plot.title = element_text(colour="black", size=panel.lab.sz, face="bold", hjust=0.5),
legend.position="none")
Tol_PAM.ES.plot<-Tol_PAM.ES.plot + ggtitle("Photo. Effic.") + labs(x="") +
theme(legend.position=c(.8, .8))
##Symbiont Retention
TP1_TolSym.plot<- TP1_TolSym.plot + labs(x="") + theme(legend.position="none")
TP2_TolSym.plot<- TP2_TolSym.plot + labs(x="", y= "") + theme(legend.position="none")
TP3_TolSym.plot<- TP3_TolSym.plot + labs(x="", y= "") + theme(legend.position="none")
TP4_TolSym.plot<- TP4_TolSym.plot + labs(x="", y= "") + theme(legend.position="none")
Tol_Sym.ES.plot<-Tol_Sym.ES.plot + ggtitle("Symbionts") + labs(x="") + theme(legend.position="none")
##Chlorophyll Retention
TP1_TolChl.plot<-TP1_TolChl.plot + theme(legend.position="none")
TP2_TolChl.plot<-TP2_TolChl.plot + labs(y="") + theme(legend.position="none")
TP3_TolChl.plot<-TP3_TolChl.plot + labs(y="") + theme(legend.position="none")
TP4_TolChl.plot<-TP4_TolChl.plot + labs(y="") + theme(legend.position="none")
Tol_Chl.ES.plot<-Tol_Chl.ES.plot + ggtitle("Chlorophyll")+ theme(legend.position="none")
##Create Panel
Tolerance_fig<-plot_grid(TP1_TolPAM.plot, TP2_TolPAM.plot,
TP3_TolPAM.plot, TP4_TolPAM.plot, Tol_PAM.ES.plot,
TP1_TolSym.plot, TP2_TolSym.plot,
TP3_TolSym.plot, TP4_TolSym.plot, Tol_Sym.ES.plot,
TP1_TolChl.plot, TP2_TolChl.plot,
TP3_TolChl.plot, TP4_TolChl.plot, Tol_Chl.ES.plot,
rel_widths=c(0.8, 0.8, 0.8, 0.8, 1), rel_heights = 1,
nrow=3, ncol=5, byrow=T, labels = NULL)
Warning: Removed 5 rows containing missing values (`geom_text()`).Warning: Removed 7 rows containing missing values (`geom_text()`).Warning: Removed 8 rows containing missing values (`geom_text()`).
##Save Figure
ggsave(filename="Figures/03_Performance/Fig5_Tolerance.png", plot=Tolerance_fig, dpi=300, width=16, height=12, units="in")
##Timepoint 1 to 2
T1.2_TLE_PAM.plot<-T1.2_TLE_PAM.plot + ggtitle("Photo. Effic.")+ labs(x="", y = "Tolerance")+
theme(plot.title = element_text(colour="black", size=panel.lab.sz, hjust=0.5), legend.position="none")
T1.2_TLE_Sym.plot<- T1.2_TLE_Sym.plot + ggtitle("Symbionts")+ labs(x="", y = "")+
theme(plot.title = element_text(colour="black", size=panel.lab.sz, hjust=0.5), legend.position="none") + ylim(-.15, .8)
T1.2_TLE_Chl.plot<- T1.2_TLE_Chl.plot + ggtitle("Chlorophyll")+ labs(x="", y = "")+
theme(plot.title = element_text(colour="black", size=panel.lab.sz, hjust=0.5), legend.position="right")
##Timepoint 3 to 4
T3.4_TLE_PAM.plot<-T3.4_TLE_PAM.plot + labs(x="Growth", y = "Tolerance") +
theme(legend.position="none")
T3.4_TLE_Sym.plot<-T3.4_TLE_Sym.plot + labs(x="Growth", y = "") +
theme(legend.position="none")
T3.4_TLE_Chl.plot<-T3.4_TLE_Chl.plot + labs(x="Growth", y = "") +
theme(legend.position="right")
##Create Panel
TradeOff_fig<-plot_grid(T1.2_TLE_PAM.plot, T1.2_TLE_Sym.plot, T1.2_TLE_Chl.plot,
T3.4_TLE_PAM.plot, T3.4_TLE_Sym.plot, T3.4_TLE_Chl.plot,
rel_widths=c(.7, .7, 1, .7, .7, 1), rel_heights = c(1, .9),
nrow=2, ncol=3, byrow=T,
labels = c("A", "B", "C", "D", "E", "F"),
label_size=panel.lab.sz)
`geom_smooth()` using formula = 'y ~ x'Warning: is.na() applied to non-(list or vector) of type 'expression'`geom_smooth()` using formula = 'y ~ x'Warning: Removed 1 rows containing non-finite values (`stat_smooth()`).Warning: Removed 1 rows containing missing values (`geom_point()`).Warning: is.na() applied to non-(list or vector) of type 'expression'`geom_smooth()` using formula = 'y ~ x'Warning: is.na() applied to non-(list or vector) of type 'expression'`geom_smooth()` using formula = 'y ~ x'Warning: is.na() applied to non-(list or vector) of type 'expression'`geom_smooth()` using formula = 'y ~ x'Warning: is.na() applied to non-(list or vector) of type 'expression'`geom_smooth()` using formula = 'y ~ x'Warning: Removed 1 rows containing non-finite values (`stat_smooth()`).Warning: Removed 1 rows containing missing values (`geom_point()`).Warning: is.na() applied to non-(list or vector) of type 'expression'
##Save Figure
ggsave(filename="Figures/03_Performance/Fig6_TradeOff.png", plot=TradeOff_fig, dpi=300, width=12, height=8, units="in")
##FvFm
IN_PAM.plot<-IN_PAM.plot + labs(x="", y="Fv/Fm") + theme(legend.position="none") + ylim(0.5, 0.655)
IN_TolPAM.plot<-IN_TolPAM.plot + labs(x="", y="Retention") +
theme(legend.position="none") + ylim(0.75, 1)
##Symbionts
IN_Sym.plot<-IN_Sym.plot + labs(x="", y="Symbionts") + theme(legend.position="none") + ylim(0,1)
IN_TolSym.plot<-IN_TolSym.plot + labs(x="", y="Retention") +
theme(legend.position="none") + ylim(.25, .75)
##Chlorophyll
IN_Chl.plot<-IN_Chl.plot + labs(x="Site and Treatment", y="Chlorophyll") +
theme(legend.position="none") + ylim(0,1.25)
IN_TolChl.plot<-IN_TolChl.plot + labs(x="Site", y="Retention") +
theme(legend.position="none") + ylim(0, 0.25)
##Create Panel
IN_HeatAssay_fig<-plot_grid(IN_PAM.plot, IN_TolPAM.plot,
IN_Sym.plot, IN_TolSym.plot,
IN_Chl.plot, IN_TolChl.plot,
rel_widths=c(1, 0.75),
rel_heights = c(1, 1, 1),
nrow=3, ncol=2, byrow=T,
labels = c("A", "B", "C", "D", "E", "F"),
label_size=panel.lab.sz)
##Save Figure
ggsave(filename="Figures/03_Performance/FigS4_InitialHeatAssay.png", plot=IN_HeatAssay_fig, dpi=300, width=8, height=11, units="in")
##Combine Results Tables
TableS5A_Growth<-data.frame(rbind(TLE_T1.2.res, TLE_T3.4.res))
TableS5B_Tolerance<-data.frame(rbind(Tol_T1.res, Tol_T2.res, Tol_T2.res, Tol_T4.res))
names(TableS5B_Tolerance)[9]<-"Timepoint"
TableS5_Perf.LM<-rbind(TableS5A_Growth, TableS5B_Tolerance)
##Organize
names(TableS5_Perf.LM)
[1] "Estimate" "Std..Error" "df" "t.value" "Pr...t.." "Predictor"
[7] "EtaSq" "Response" "Timepoint"
TableS5_Perf.LM<-TableS5_Perf.LM %>% dplyr::rename(SE = Std..Error, DF = df, t = t.value, p = Pr...t.., EffectSize = EtaSq)
TableS5_Perf.LM<-TableS5_Perf.LM[,c("Timepoint", "Response", "Predictor", "Estimate", "SE", "DF", "t", "p", "EffectSize")]
#Round to 3 digits
TableS5_Perf.LM$Estimate<-round(TableS5_Perf.LM$Estimate, 3)
TableS5_Perf.LM$SE<-round(TableS5_Perf.LM$SE, 3)
TableS5_Perf.LM$t<-round(TableS5_Perf.LM$t, 3)
TableS5_Perf.LM$p<-round(TableS5_Perf.LM$p, 3)
TableS5_Perf.LM$EffectSize<-round(TableS5_Perf.LM$EffectSize, 3)
#Integer
TableS5_Perf.LM$DF<-round(TableS5_Perf.LM$DF, 0)
##Write Out Table
write.csv(TableS5_Perf.LM, "Tables/TableS5_Performance_LM_Results.csv", row.names=FALSE)
##Combine Results Tables
TableS6_IN_Tolerance<-data.frame(rbind(PAM_IN.res, Sym_IN.res, Chl_IN.res))
##Organize
names(TableS6_IN_Tolerance)
[1] "estimate" "SE" "df" "t.ratio" "p.value" "Predictor"
[7] "Response" "EtaSq"
TableS6_IN_Tolerance<-TableS6_IN_Tolerance %>% dplyr::rename(Estimate = estimate, DF = df, t = t.ratio, p = p.value, EffectSize = EtaSq)
TableS6_IN_Tolerance<-TableS6_IN_Tolerance[,c("Response", "Predictor", "Estimate", "SE", "DF", "t", "p", "EffectSize")]
#Round to 3 digits
TableS6_IN_Tolerance$Estimate<-round(TableS6_IN_Tolerance$Estimate, 3)
TableS6_IN_Tolerance$SE<-round(TableS6_IN_Tolerance$SE, 3)
TableS6_IN_Tolerance$t<-round(TableS6_IN_Tolerance$t, 3)
TableS6_IN_Tolerance$p<-round(TableS6_IN_Tolerance$p, 3)
TableS6_IN_Tolerance$EffectSize<-round(TableS6_IN_Tolerance$EffectSize, 3)
#Integer
TableS6_IN_Tolerance$DF<-round(TableS6_IN_Tolerance$DF, 0)
##Write Out Table
write.csv(TableS6_IN_Tolerance, "Tables/TableS6_Initial_Tolerance_LM_Results.csv", row.names=FALSE)